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ABSTRACT 


§  This  report  addresses  the  design  of  a  spread  spectrum  image  transmission  system 

to  provide  increased  antijam  protection  to  a  television  link  from  a  small  remotely  piloted 
I'  vehicle.  Previous  quarterly  reports  have  described  component  developments,  theoretical 

advances  in  image  processing,  and  simulations  of  proposed  coding  schemes.  This  report 
J  summarizes  the  proposed  total  system. 
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REPORT  SUMMARY 


This  report  describes  the  proposed  system  for  bandwidth  reduction  o*f  the  video 
for  a  remotely  piloted  vehicle.  Two  systems  ?re  being  constructed:  one  compatible  with 
a  standard  vidicon,  the  other  compatible  with  a  Charge  Injection  Device  100  x  100  camera. 
The  transversal  filters  for  the  100  point  Cosine  Transform  were  delivered  during  this  phase 
and  integrated  into  a  transform  subsystem.  A  limited  prototype  ground  station  was  de¬ 
signed  and  constructed  during  this  phase.  The  proposed  transform  was  shown  to  be  rela¬ 
tively  tolerant  to  bit  errors,  which  is  shown  in  Appendix  A.  The  proposed  coding  schemes 
were  investigated  and  are  reported  in  Appendix  B.  The  modularity  of  the  CHIRP-Z  Trans¬ 
form  is  shown  in  Appendix  C.  This  modularity  and  the  structures  described  in  Apper.dix  D 
allow  general  purpose  high-speed  signal  processors  which  could  be  used  in  the  channel 
decoder. 
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INTRODUCTION 


The  Image  Transmission  via  Spread  Spectrum  Techniques  program  has  developed 
to  a  point  where  a  complete  proposed  system  can  be  described.  Previous  quarterly  reports 
have  described  component  developments,  theoretical  advances  in  image  processing,  and 
simulations  of  proposed  coding  schemes.  A  summary  of  the  complete  system  description 
was  presented  during  an  ARPA  workshop  at  MITRE  on  May  14-15,  1974.  This  summary 
is  presented  below.  Also  included  are  Appendixes  A  thrugh  F  which  present  study  efforts 
conducted  during  this  last  reporting  period  by  individual  authors. 

Appendix  A  is  a  compilation  of  photographs  of  transformed  images  done  at  the 
University  of  Southern  California.  Appendix  B  is  a  comparison  of  methods  of  channel 
encoding.  Appendix  C  is  a  description  of  the  modular  CH1RP-Z  Transform.  Appendix  D  is 
a  description  of  a  high  capacity  linear  processor.  Appendix  E  is  a  description  of  the  utiliza¬ 
tion  of  a  signal  processing  imager  chip.  Appendix  F  is  an  updated  description  of  a  signal 
processing  interpreter  for  use  in  a  general  purpose  computer. 


SYSTEM  DESCRIPTION 


The  following  material  is  presented  in  summary  format  and  represents  the  status 
of  the  complete  proposed  program. 

Problem  Statement 

Design  and  demonstrate  a  real  time  image  bandwidth  reduction  and  spread  spectrum 
encoding  system  for  use  in  a  lightweight  Remotely  Piloted  Vehicle  (RPV)  television  system 
to  provide  protection  against  intentional  and  unintentional  jamming. 

Goals  for  RPV  Image  Redundancy  Reduction  System 

a)  -  1.51b 

b)  15 -30  watts 

c)  Field  demonstration  in  spring  1975. 

Approach 

Reduce  data  rate  in  order  to  use  spread  spectrum  encoding 

•  Number  of  image  points 

•  Frame  rate 

•  Transform  encoding 

Performance  of  Intraframe  Encoding 

PCM/DPCM 
8-3  bits/pixel 
Transform  with 

Selective:  3  -  1  bits/pixel 
Quantization 
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Transforms  Examined: 


•  Haar 

•  Hadamard 

•  Slant 

•  Fourier 

•  Differential  Pulse  Code  Modulation  (DPCM) 

•  Discrete  Cosine  Transform  (DCT) 

•  Karhunen  -  Loeve  (K  -  L) 


Mean  Square  Error  Performance  of  Various  Transforms  for  Exponentially 
Correlated  Data,  p  =  0.9: 


Discrete  Fourier  Transform,  CZT  Algorithm: 


HG  =  IT1  e-j27rmn/N  „ 

U  m  nV  8n 


2)  -2mn  =  -rrr  +  (m  -  n)“  -  n“ 


3)  :.G 


m 


where  m  =  0,  1,  2, .  . . ,  N  -  1. 


CHIRP-Z  Transform  Implementation  of  the  DFT 


★  Denotes  either  convolution  or  circular  convolution 
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DFT  Via  CZT  Algorithm  with  Parallel  Implementation  of  Complex  Arithmetic 


Relation  of  DCT  to  DFT 


DCT  -  DFT  of  symmetric  extension  =  real  part  of  DFT  of  null  extension 


Serial  Access  Implementation  of  the  DCT 


.Jirn'/ON-I) 


Charge  Transfer  Implementation  of  DCT 


CHARGE  TRANSFER 
_ LSI  CHIP 


sin  7rmJ/(2N-1) 


Transforms  Examined 


DPCM  Encoder 


Two  -  Dimensional  Hybrid 


Haar 

Haar/DCT 

Hadamard 

Hadamard/DPCM 

Fourier 

Fourier/DPCM 

DPCM 

DCT/DPCM 

DCT 

MEMORY 

CONTROL 

STATE 

MACHINE 

10s  Samples/Second  Transform  System 


•  100  X  100  CID  sensor 

•  100  point  BBD  implementation  of  DCT/DPCM 


5  X  106  Samples/Second  Transform  System 

•  256  X  240  sensor 

•  32  point  CCD  implementation  of  DCT/DPCM 


Coding  Design 

Source 

a)  minimum  of  4  frames/sec 

b)  256  X  256  pixels 

c)  mean  square  error  *  30  dB  down 
Channel 

a)  20  Mbits/sec 

b)  at  each  data  rate  maximize  distance  between  code  words 

c)  coding  for  modem  independence 


System  Performance  for  a  Linear  Matched  Filter  PN-PPM  Receiver 


FRAMES/ 


SEC 

BITS/PIXEL 

30 

1 

.56 

.31 

15 

1.75 

1 

.56 

.31 

7-1/2 

3 

1.75 

1 

.56 

.31 

3-3/4 

3 

1.75 

1 

.56 

1-7/8 

3 

1.75 

1 

15/16 

3 

1.75 

18 

21 

24 

27 

30 

MATCHED  FILTER  PROCESSING  GAIN 


Performance  of  Proposed  System 


GAUSSIAN  NOISE  PULSED 


CODE 

WORD 

BITS 

PER 

COMPRESSION 

RATIO 

MARGIN 

INTERFERENCE 

MARGIN 

LENGTH 

PIXEL 

P=  10'J 

P  =  10'3 

(CODEWORD 

SEPARATION) 

1023 

0.31 

102.3 

22.4  dB 

23.6  dB 

30.1  dB 

511 

0.56 

56.8 

19.7  dB 

20.5  dB 

27.1  dB 

255 

1.0 

32.0 

17.2  dB 

18.1  dB 

24.1  dB 

*P=  10's 
S/N  *  10  dB 


GAUSSIAN 

NOISE 

MARGIN 

OVER  SPEC* 


27.7  dB 
25.0  dB 
22.5  dB 


Transmitter 


SOURCE  ENCODER  CHANNEL  ENCODER 

f - - N  ( - ~ N 


Receiver 


IF 


CHANNEL  DECODER 


SOURCE  DECDDER 


DISPLAY 


APPENDIX  A 


UNIVERSITY  OF  SOUTHERN  CALIFORNIA 
TRANSFORM  SIMULATION 


TRANSFORMED  IMAGES 


This  appendix  presents  a  compilation  of  photographs  of  several  image  coding  systems 
based  upon  transform  coding  performed  at  the  University  of  Southern  California. 
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2  Dim.  Cosiru*, 


Cosine  /  DPv^M, 


Cosin e/DFCM,  P  -  10 


2  Dim.  Cosine,  P  =  10 


Original  Hadamard  /  DPCM,  P  0 

Figure  1.  Coded  Pictures,  2  Bits/Pixel. 


Cosine/ DPCM,  P  0 


2  Dim.  Cosine,  P  0 


Cosine/DPCM,  P  10 


2  Dim.  Cosine,  P  =  10 


Hadamard/DPCM,  P  -  0 


Figure  2.  Coded  Pictures,  1  Bit/Pixel. 


Cos/DPCM,  P 


Cos/DPCM,  P 


Cos /D  PCM,  P 


Cos/DPCM 


Cos/DPCM,  P  10 


Cos/DPCM,  P 


Figuro  6.  Coded  Pictures,  1  Hit/Pixel 


Cos/DPCM 


Cos/DPCM 


Cos/DPCM,  P  =  10 


Cos/DPCM,  P 


Figure  7.  Coded  Pictures,  0.5  Bit/Pixel 


original 


Cosim'/DPCM,  P 


Cosi nc / I)PCM ,  P  10"  Cosine/DPCM,  P 


Figure  9.  Coded  Pictures,  2  Bits/Pixel. 


* 


Cos/DPCM,  P 


Cos/DPCM,  P 


Cos/DPCM,  P  =  10 


Cos/DPCM,  P 


Figure  10.  Coded  Pictures,  1  Bit/Pixel 


*»% 

4  < 

r*  7J 

4k. 

■BF 

n 

•<  J 

Cos/DPCM,  P  =  0 


Cos/DPCM,  P  0 


Cos/DPCM,  P  =  10 


Cos/DPCM,  P  =  10 


Cos/DPCM,  P  =  10 


Cos/DPCM,  P  =  10 
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VIDEO  SPREAD  SPECTRUM  ENCODING 
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Dr.  E.  H.  Wrench,  Jr. 

Naval  Undersea  Center 
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VIDEO  SPREAD  SPECTRUM  ENCODING 


by 

Edwin  H.  Wrench,  Jr. 


INTRODUCTION 


In  the  design  of  a  communication  system  there  is  usually  a  trade  off  between  system 
performance  and  complexity  (expense).  Such  trade  offs  between  complexity  and  perform¬ 
ance  in  the  presence  of  noise  exist  for  the  video  down  link  for  the  ARPA  remotely  piloted 
vehicle.  The  allocated  channel  bandwidth  is  barely  adequate  to  transmit  the  video  directly; 
however,  if  direct  transmission  is  used,  noise  immunity  can  only  be  obtained  by  increased 
transmitted  power.  I'or  a  system  of  limited  si/e.  weight,  and  power  dissipation,  this  approach 
is  undesirable. 

An  alternative  approach  is  to  use  bandwidth  compression  techniques  to  remove  re¬ 
dundancy  from  the  video  prior  to  transmission.  The  compressed  video  can  then  be  trans¬ 
mitted  directly  or  recoded  to  fill  the  allocated  channel  bandwidth.  The  latter  technique 
further  increases  noise  immunity  by  trading  bandwidth  for  increased  signal-to-noise  ratio. 
This  more  complicated  system  offers  the  same  performance  as  the  direct  transmission  scheme 
but  with  reduced  transmitted  power. 

In  this  paper  the  theoretical  performance  of  three  proposed  video  down  link  com¬ 
munication  systems  are  compared  with  the  performance  of  a  reference  system  which  uses 
direct  transmission  of  the  raw  video  data.  The  first  system  (system  1 )  transmits  the  binary 
compressed  video  directly  using  biorthogonal  code  words.  The  second  system  (system  2) 
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uses  M-ary  signaling  techniques  to  increase  the  bandwidth  of  the  compressed  video  to  that 
of  the  original  raw  video  data.  The  code  words  are  pseudorandom  noise  (PN)  sequences 
generated  by  maximal  length  shift  register  code  generating  techniques.  They  have  the  prop¬ 
erty  that  the  auto  correlation  function  has  constant  side  lobes  of  -1/M,  where  M  is  the 
length  of  the  code  word.  In  the  receiver,  the  received  signal  is  correlated  with  M  code  refer¬ 
ences.  The  correlator  with  the  largest  output  is  chosen  as  the  transmitted  code  word 
(matched  filter  detection).  The  third  system  (system  3)  uses  code  words  that  are  identical 
to  those  of  system  2,  but  uses  binary  transmission  with  biorthogonal  code  words.  The  M 
bits  of  the  code  words  are  then  detected  on  a  bit-to-bit  basis  and  the  resulting  detected  bit 
stream  fed  into  the  correlators.  As  with  system  2,  the  code  word  corresponding  to  the 
largest  correlator  output  is  chosen  as  the  transmitted  code  word.  The  difference  between 
systems  2  and  3  is  that  system  3  detection  is  performed  prior  to  correlation  rather  than 
using  correlation  to  perform  the  detection  as  in  system  2. 

The  performance  of  all  three  systems  will  be  compared  with  that  of  the  reference 

system  (system  4),  which  transmits  and  detects  a  20  megabit  per  second  video  on  a  bit-to- 

bit  basis.  The  code  words  used  are  biorthogonal.  The  noise  is  assumed  to  be  Gaussian  dis- 

-2 

tributed  white  noise  in  all  cases  and  the  bit  error  rate  taken  as  10  . 

Several  different  compression  ratios  for  the  video  were  examined.  The  compression 
ratios  for  various  length  PN  sequences  and  for  a  frame  rate  reduction  factor  of  eight  are 
given  in  Table  1. 

Table  1.  Compression  Ratio  for  Various  Code  Lengths 
for  a  Frame  Rate  Reduction  Factor  of  Eight 


PN  Length,  " M  " 

Compression  Ratio 

255 

32 

511 

56.8 

1023 

102.4 

Analysis  of  the  Reference  System  (System  4) 


The  signal-to-noise  ratio  7  at  the  output  of  the  matched  filter  required  to  achieve 
a  given  bit  error  rate  in  a  binary  system  with  biorthogonal  code  words  is 

7=  (erfc-1  (2Pb)]2. 


For  a  bit  error  rate  (P^)  of  10 


7  «  2.64  =  4.3  dB. 


Analysis  of  System  1 

System  1  is  identical  to  the  reference  system  (system  4)  except  the  transmission 
rate  has  been  reduced  by  a  compression  ratio  factor  (C),  where 


_  No.  of  bits/sec  original  data 
~  No.  of  bits/sec  compressed  data 


The  output  signal-to-noise  ratio  of  the  matched  filter  is 


r  setidt 

JL  N0  ’ 


where  S(t)  is  the  instantaneous  received  signal  power,  and  N0  is  the  noise  spectral  density. 
Since  the  duration  of  the  signal  pulse  has  increased  for  system  1  by  a  factor  of  C, 


f  S.(t)dt  =  cf 

*/— 00  *  J-oc 


S4<t)dt. 
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Therefore  the  signal-to-noise  ratio  for  system  1  is  increased  by  a  factor  of  C  over  that  of 
system  4.  The  improvement  in  signal-to-noise  ratio  for  3  values  of  C  are  given  in  Table  2. 


Table  2.  Signal-to-Noise  Ratio  Improvement  (dB)  for  System  1 


Compression 
Ratio  f  C) 

Improvement 

(dB) 

32 

15.1 

56.8 

17.5 

102.4 

20.1 

Analysis  of  System  2 


The  performance  of  System  2,  employing  M-ary  signaling  with  matched  filter  de¬ 
tection,  is  described  by  the  probability  of  a  symbol  error, 


QQ  2 

P  =  1  -  f  (2i T)'1/2  e"x  /2 


x+(2t) 


1/2 

(27T)-1/2 


m-1 


dx. 


For  large  M,  probability  of  a  bit  error  is 


P 


b 


=  1/2  P  . 
s 


The  values  of  P  for  various  values  of  M  and  y  are  tabulated  by  Viterbi.  * 

The  improvement  in  the  required  signal-to-noise  ratio  is  given  in  Table  3  for  various 
compression  ratios. 

1  Viterbi,  Andrew  J.  Phase-Coherent  Communication  over  Continuous  Gaussian  Channel  in  Digital  Communications  with 
Space  Applications  by  Solomon  W.  Golomb.  Prentice  Mall,  1964.  ~ 
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Therefore  the  signal-to-noise  ratio  for  system  1  is  increased  by  a  factor  of  C  over  that  of 
system  4.  The  improvement  in  signal-to-ioise  ratio  for  3  values  of  C  are  given  in  Table  2. 


Table  2.  Signal-to-Noise  Ratio  Improvement  (dB)  for  System  1 


Compression 
Ratio  ( C) 

Improvement 

(dB) 

32 

15.1 

56.8 

17.5 

102.4 

20.1 

Analysis  of  System  2 

The  performance  of  System  2,  employing  M-ary  signaling  with  matched  filter  de¬ 
tection,  is  described  by  the  probability  of  a  symbol  error, 

f  \rH2y)'nw-mt-y2Riy 

L«/_oo 

For  large  M,  probability  of  a  bit  error  is 


The  values  of  Ps  for  various  values  of  M  and  7  are  tabulated  by  Viterbi. 1 

The  improvement  in  the  required  signal-to-noise  ratio  is  given  in  Table  3  for  various 
compression  ratios. 

*  Viterbi,  Andrew  J.  Phase-Coherent  Communication  over  Continuous  Gaussian  Channel  in  Digital  Communications  with 
Space  Applications  by  Solomon  W.  Goiomb.  Prentice  Hall,  1964.  ~  — 


m-1 

dx. 
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Table  3.  Signal-to-Noise  Ratio  Improvement  (dB)  for  System  2 


Compression 

Improvement 

Ratio  (C) 

(uB) 

32 

17.2 

56.8 

19.7 

102.4 

22.4 

Analysis  of  System  3 

System  3  uses  matched  filtering  following  detection  on  a  bit-to-bit  basis.  The  prob¬ 
ability  of  a  symbol  error  is  derived  as  follows:  There  are  M  correlators,  one  for  each  code 
word.  The  output  from  the  correlator  which  corresponds  to  the  transmitted  code  has  a 
value  of 


=  M  -  2E 

when  exactly  E  bit  errors  have  occurred.  The  correlators  with  only  side  lobe  sums  have 
values  of 


SK  =  -I  -2A  +  2D, 

where  A  is  the  number  of  bit  errors  in  bits  that  should  agree  with  the  reference,  and  D  is 
the  number  of  bit  errors  in  bits  that  should  disagree  with  the  reference.  Since  A  +  D  =  E, 
the  equation  can  be  rewritten 


SL  =  -1  +  4D-2E. 
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For  an  individual  correlator,  an  error  will  occur  if  the  side  lobe  sum  exceeds  the  main  lobe 
sum,  i.e.,  S£  >  Mg.  There  is  also  a  possible  error  when  S£  =  M£.  The  probability  of 
making  an  error  in  any  correlator,  given  E  bit  errors,  is 


Pc/E  =  Prob  (SE  >  Me)  +  1/2  Prob  (S£  =  M£). 


M+l. 


Cut,  Prob  (SP  >  Mp)  =  Prob  (-1-4A-2E  >  M-2E)  =  Prob  (D  > 


No.  errors  =  E 


so  that 


M+k 


M+l. 


Pc/E  =  Prob  (D  >  ^)  £  +  1/2  Prob  (D  = 


E’ 


or,  p  =l-HYGEO(^,M,^,E)+l/2[HYGEO(^,M,^-,E) 


M+l  „  M-l 


c/E 


-HYGEO  (^pM,^,  E», 


where  HYGEO  is  the  hypergeometric  distribution 

H«K.N,N,,NR)=g(^)f"Ry("i). 

If  we  assume  the  probability  of  an  error  in  a  correlator  is  independent  of  the  error  in  the 
other  correlators,  the  total  probability  of  a  correlator  error  is 

P  =  (1  -Prob.  of  no  correlator  errors)  *  1- (1  ’Pc/E^^  ^ 
c 


The  probability  of  a  symbol  error  can  then  be  written  as 

M  M 

Ps=  V  Prob  (E)*P  =  Yj  Pr°b(E)*Pc, 

E=0  E=[(M+l)/4) 

where  Prob  (E)  is  the  probability  of  E  errors 


Prob  (E)  =  (P)E(1-P)M'E(e), 
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and  where  P  is  the  probability  of  a  single  bit  error  (determined  by  S/N)  in  the  detector. 
For  biorthogonal  coherent  detection 

P  =  1/2  crfc  \fy  . 

Solution  of  the  equation  for  Ps  yields  the  results  shown  in  Figure  1.  For  M  =  255,  and 

P,  =  0.01  =  1/2  P. ,  the  required  signal-to-noise  ratio  is  -10.8  dB.  This  same  bit  error  prob- 
b  b 

ability  is  obtained  without  bandwidth  compression  and  coding  with  S/N  *  4.2  dB.  The 
improvement  is  then  15.0  dB  for  the  compressed  coded  video  over  straight  video  trans¬ 
mission.  These  results  and  those  for  M  =  5 1 1  and  M  =  1023  are  given  in  Table  4. 

_2 

Table  4.  Signal-to-Noise  Ratio  Improvement  (dB)  for  System  3  < Pj  =10*") 


Compression 
Ratio  ( C) 

Improvement 

(dB) 

32 

15.5 

56.8 

18.0 

102.4 

20.8 

Discussion 

The  analysis  results  are  summarized  for  the  three  systems  in  Table  5  for  a  probability 
of  a  bit  error  equal  to  10  *\  Similar  results  are  given  for  a  probability  of  a  bit  error  equal 
to  10'3  in  Table  6.  Table  7  summarizes  the  improvement  for  System  2  with  pulse  inter¬ 
ference  rather  than  Gaussian  noise. 

Each  of  the  three  systems  examined  has  advantages  and  disadvantages.  System  1 
gives  a  15-20  dB  improvement  over  the  reference  system  (system  4)  and  requires  no  sophis¬ 
ticated  coding  or  decoding.  However,  it  is  a  relatively  narrow  band  system,  depending  on 
the  compression  ratio,  that  lends  itself  to  easy  detection  by  a  jammer. 
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System  2  gives  a  17.2  -  23.7  dB  improvement,  and  wide  bandwidth,  at  the  expense 
of  greatly  increased  system  complexity,  since  it  requires  matched  filters  at  IF  frequencies. 
System  3  is  a  compromise  between  system  1  and  system  4.  The  transmission  and  detection 
are  done  one  bit  at  a  time  to  eliminate  the  need  for  large  analog  matched  filters.  The  per¬ 
formance  is  within  2  dB  of  the  optimum,  provided  by  system  2. 

Table  5.  Summary  of  System  Improvement  (dB)  with  the  Probability 

.2 

of  a  Bn  Error  Equal  to  10  and  Gaussian  Noise 


System 

Compression  Ratio  (C) 

32  56.8  102.4 

1 

15.1 

17.5 

20.1 

2 

17.2 

19.7 

22.4 

3 

15.5 

18.0' 

20.8 

Table  6.  Summary  of  System  Improvement  (dB)  with  the  Probability 

-3 

of  a  Bit  Error  Equal  to  10  and  Gaussian  Noise 


System 

Compression  Ratio  ( C) 

32  56.8  102.4 

1 

15.1 

17.5 

20.1 

2 

18.4 

21.0 

23.7 

3 

16.4 

19.3 

22.1 

Table  7.  System  2  Improvement  (dB)  with  Pulse  Noise 
(10  Log  [Code  Word  Separation) ) 


Compression  Ratio  ( C) 

32 

56.8 

102.4 

Improvement  (dB) 

21 

24 

27 
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Probability  of  a  Code  Word  Error  (Ps) 


APPENDIX  C 

THE  MODULAR  CHIRP-Z  TRANSFORM 

J.  M.  Speiser  and  H.  J.  Whitehouse 

Naval  Undersea  Center 


ABSTRACT 


Two  methods  for  combining  P  Chirp-Z  Transform  (CZT)  modules  of  length  N  to 
perform  a  Discrete  Fourier  Transform  (DFT)  of  length  NP  are  described.  The  firs 
method  uses  an  auxiliary  parallel-input,  parallel-output  DFT  device  of  size  P 
and  allows  the  transform  of  size  NP  to  be  performed  in  the  same  time  required  fc 
a  single  CZT  module  to  perform  a  size  N  transform.  The  second  method  uses  an 
auxiliary  parallel -input,  serial-output  DFT  device  of  size  P.  If  the  second 
method  is  implemented  entirely  in  a  single  technology,  such  as  with  CCD's,  it 
performs  the  size  NP  transform  in  P  times  the  amount  of  time  required  for  a 
single  CZT  module  to  perform  a  size  N  transform;  if  P  is  a  composite  number, 
say  P  =  Pi P2 »  the  second  method  also  permits  the  same  hardware  to  perform  P^ 
simultaneous  transforms  of  length  NP2- 
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INTRODUCTION 


Many  signal  processing  problems  require  flexible  real-time  implementations 
of  linear  signal  processing  operations  such  as  Fourier  transforms,  convolution 
correlation,  and  beamforming.  All  of  these  operations  may  be  performed  at  high 
throughput  rates  using  the  discrete  Fourier  transform  implemented  via  the  CZT 
algorithm  with  a  transversal  filter  or  cross  convolver  used  to  perform  the 
required  convolution  or  correlation  with  a  discrete  chirp  [1-3]  as  shown  in 
Fig.  1  and  described  by  equations  (1)  -  (3).  Special  purpose  methods  have  also 
been  previously  described  for  combining  a  number  of  CZT  modules  to  perform  a 
longer  DFT  [4].  The  previous  methods,  however,  required  a  different  acoustic 
surface  wave  filter  for  each  different  number  of  CZT  modules  to  be  combined. 

This  limits  the  flexibility  of  the  transform  size  attainable  with  a  given  set 
of  components,  and  also  prevents  the  longer  transform  system  from  being  externally 
clocked,  since  the  propagation  time  through  the  surface  wave  device  cannot  be 
varied  by  more  than  a  small  fraction  of  one  percent. 

A  one-dimensional  discrete  Fourier  transform  may  be  written  as  a  partial 
transform  of  a  doubly  subscripted  representation  of  the  data,  followed  by  a 
pointwise  multiplication,  followed  by  a  second  partial  transform  [4,5]  as  shown 
in  equations  (4)  -  (9). 

M-1 

Hk  ■  £hi 

j=0 


-i27rjk 

M 


(1) 
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(2) 


M-1 


H, 


pk  2  ihi  pri  pk-i 

i=o 


ps  =  e 


I7TS' 

M 


=  P 
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pn-1  -iMn 

PN 
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Ee 
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-i2Tik(Pn+s) 

PN 


PN 


£ 
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9Pn+s 


-i2jrkji 

N 


Let  Q$(k)  =  e 


-J2irJkS 
PN 


Q$(k+nN)  =  e 


-i2Ti(k+nN)s 
PN 
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a  PN  PN 

e  e 


PN-1 
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(5) 


(6) 


(7) 
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-i  2iTsn 


Q$(k+nN)  =  e 


Qs00  =  csn  qs(k) 


-i  2Tisn 


where  c  =  e 
sn 


(9) 


Figures  2  and  2A  show  modular  CZT  implementations  which  follow  from 
equations  (4)  -  (9).  The  individual  CZT  subsystems  shown  in  Figures  2  and  2A 
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would  be  similar  to  the  CZT  implementations  previously  described  [1-3].  A 
suitable  serial  to  parallel  multiplexer  is  shown  in  Fig.  3. 

The  parallel  DFT  required  for  the  second  partial  transform  in  Fig.  2  may  be 
implemented  as  combination  of  summers  and  attenuators.  This  is  shown  in  Figures 
4-6  for  P  =  2,  3,  and  4.  In  general,  the  attenuation  factors  are  complex, 
requiring  separate  weightings  for  the  real  and  imaginary  parts  as  shown  in  Fig. 
7.  A  complete  double  length  CZT  is  shown  in  Fig.  8.  Unfortunately,  a  parallel 
DFT  implementation  of  this  type  becomes  unwieldy  if  the  dimension  P  is  very 
large. 

Other  types  of  discrete  Fourier  transform  implementations  may  be  derived 
from  the  identities  of  equations  (10)  -  (12). 

P-1  -i_2jrkn 

Hk  *  S e  p  hk 


kn  =  J-((k+n)2  -  (k-n)2) 


These  equations  have  been  used  by  Means  [6]  to  design  a  transform  device  in  which 
signals  are  shifted  through  two  delay  lines  at  different  speeds.  Alternatively, 
if  the  factors  in  equation  (12)  are  interpreted  as  two  waves  propagating  in 
opposite  directions  relative  to  the  function  to  be  transformed,  it  may  be  seen 
that  the  structure  of  Fig.  9  also  performs  a  discrete  Fourier  transform  with 
speed  comparable  to  that  of  a  CZT.  By  changing  the  reference  functions  applied 
to  the  delay  lines  and  partitioning  the  input  leads  and  output  summer  as  shown 


in  Fig.  10,  the  hardware  of  Fig.  9  may  also  perform  several  shorter  discrete 
Fourier  transforms.  If  the  parallel  input-serial  output  CZT  of  Fig.  10  is  used 
in  the  modular  CZT  of  Fig.  2A,  then  discrete  Fourier  transforms  of  length  P£N 
may  be  performed  simultaneously,  where  P  =  P-^  is  any  factorization  of  P^ . 
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Figure  2A.  Alternate  Organization  of  Modular  CZT 


0 


RECIRCULATING 
BINARY  SHIFT 
REGISTER 
P  BITS  LONG 


OUTPUTS 


Figure  3.  Serial  to  Parallel  Multiplexer 
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I 
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INTRODUCTION 


High  speed  serial  access  discrete  Fourier  transform  (DFT)  devices  have 
been  built  using  the  chirp-z  transform  (CZT)  with  transversal  filler  implemen¬ 
tation  of  the  required  convolution  [1-3].  Such  devices  perform  a  DFT  by 
premultiplying  by  a  discrete  chirp,  convolving  with  a  discrete  chirp,  and 
pos tmulti plying  by  a  discrete  chirp.  One  transform  point  is  obtained  per  shift, 
and  the  size  of  the  transform  is  proportional  to  the  number  of  taps  in  the 
transversal  filter  used  to  perform  the  convolution.  If  the  data  is  not 
recirculated  through  the  filter,  then  2N-1  taps  are  required  to  implement  an 
N  point  DFT.  Using  present  acoustic  surface  wave  or  CCD  transversal  filter 
technology,  it  is  difficult  to  implement  a  transversal  filter  with  more  than 
about  1000  independent  taps,  and  therefore  a  direct  transversal  filter  implemen¬ 
tation  of  the  DFT  via  a  one-dimensional  CZT  is  limited  to  a  size  of  about  500 
points . 

This  note  describes  a  DFT  architecture  based  on  a  multidimensional  CZT. 

This  architecture  uses  a  small  number  of  high  speed  transversal  filters, 
together  with  a  larger  number  of  low  speed  transversal  filters  to  implement  a 
CZT  whose  speed  is  determined  by  the  shift  rate  of  the  high  speed  transversal 
filters  and  whose  transform  size  is  determined  by  the  total  number  of  taps  of 
all  the  low  speed  transversal  filters. 
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REPRESENTATION  OF  A  ONE- DIMENSIONAL  DFT  AS  A  TWO-DIMKNS I ONAL  DFT 

It  lias  previously  been  shown  [4-6]  that  if  N  a  N^N2  w*,ere  ^  j  Ul1^  'iav0 
no  common  divisor,  then  the  one-dimensional  DFT  of  equation  (1)  and  the  two- 
dimensional  DFT  of  equation  (2)  are  equivalent. 


-12’  N1 


f  -  >  f  e 

q  U  p 


Nj-1  N2-l 


S(kl*k2)  *  zh  Xrf8(jl,J2)  6 


for  p,q  =  0,  1,  .  .  . ,  N-l 


^lkl  ^2k2 

-i2lT  n  ~  +  ~iT~ 

1  2  . 


V°  V° 


for  j^,k^  “  0,  1,  •  .  N^-l 


The  scan  required  to  make  the  two  transforms  exactly  equivalent  is  not  a 
simple  lexical  scan,  but  rather  a  linear  congruential  scan  such  as  that  shown 
in  equations  (3)  and  (4). 


”  ^lN2  +  ^2N1 


(mod  N) 


q(kltk2)  -  kiU!N2  +  k2U2Nl  (m°d  N) 


The  constants  and  U2  are  the  solutions  of  equations  (5)  and  (6) 


Vi  ‘  1 


(mod  N^) 


N1U2  “  1 


(mod  N2) 
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In  special  cases  the  output  scan  is  very  similar  to  the  input  scan.  For 
example,  if  =  N^+l,  then  U^al  and  Uj®-!,  so  that  the  two  scans  are  given 
by  equations  (7)  and  (8). 

P  =  j  +  ^2N1  (mod 

q  =  kiN2  “  k2Nl  (mod  N> 


(7) 

(8) 


.  i 
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SCAN  CONVERSION 


A  linear  congruential  scan  is  illustrated  in  figure  1,  and  a  means  of 
converting  the  one-dimensional  data  to  the  required  two-dimensional  format  is 
shown  in  figure  2.  The  multiplexer  of  figure  2  may  be  implemented  with  an 
acoustic  surface  wave  delay  line,  together  with  auxiliary  high  speed  analog 
switches.  The  most  suitable  analog  switches  at  the  present  time  would  appear 
to  be  balanced  mixers,  since  they  are  capable  of  switching  in  1-2  nanoseconds. 
As  shown  in  figures  2  and  3,  the  data  is  recirculated  once.  It  will  be  noted 
that  if  Nj^  and  N2  are  nearly  equal,  the  required  number  of  taps  is  about  N'5, 
so  that  for  a  long  transform  the  required  number  of  taps  is  much  smaller  than 
the  number  of  samples  stored  in  the  line.  This  greatly  reduces  the  problem  of 
dispersion  produced  by  acoustic  scattering  from  the  taps.  For  example,  only 
100  taps  would  be  needed  for  the  multiplexer  of  a  transform  of  size  about 
10,000.  Transversal  filters  with  100  taps  and  delay  lines  storing  10,000 
samples  are  well  within  the  existing  state  of  the  art. 
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TWO  DIMENSIONAL  CZT 


Once  the  data  is  in  a  two-dimensional  format,  with  simultaneous  serial 
access  to  all  the  rows,  it  may  be  transformed  in  the  "horizontal"  direction  by 
the  structure  shown  in  figure  4.  The  individual  one-dimensional  discrete  chirp 
filters  and  discrete  chirp  generator  of  figure  4  may  be  implemented  using  CCDs 
or  digital  correlators. 

An  acoustic  surface  wave  device  with  multiple  input  taps  may  be  used  to 
access  a  column  of  the  partially  transformed  output  in  a  single  shift  time  of 
the  partial  transform  device.  With  appropriate-  coding  of  the  surface  wave 
column  access  device,  it  may  also  perform  the  discrete  chirp  premultiplication 
and  discrete  chirp  convolution  of  a  DFT  in  the  "vertical"  direction.  A  complete 
two-dimensional  CZT  architecture  is  shown  in  figure  5,  and  the  required  coding 
for  the  column  access  surface  wave  device  is  shown  in  figure  6.  The  complex 
arithmetic  may  be  implemented  as  described  previously  [1], 

A  balanced  mixer  may  be  used  for  the  fast  multiplier  required  for  the 
vertical  transform.  Lower  speed  variable  transconductance  multipliers  may  be 
used  in  the  horizontal  partial  transform. 


HIGH  CAPACITY  DFTs 


The  scan  conversion  multiplexer  of  figure  2  may  be  combined  with  the  two- 
dimensional  CZT  of  figure  5  to  yield  a  high  capacity  CZT  device  as  shown  in 
figure  7.  The  throughput  rate  of  the  resulting  CZT  device  is  proportional  to 
the  sample  rate  of  the  surface  wave  components,  and  the  transform  size  is 
proportional  to  the  total  number  of  taps  of  all  the  transversal  filters  or 
digital  correlators  used  in  the  two-dimensional  partial  CZT. 

The  high  capacity  CZT  device  can  be  used  as  a  high  speed,  high  resolution 
spectrum  analyzer,  or  as  a  high  bandwidth  matched  filter  with  many  degrees  of 
freedom  and  variable  reference  function.  The  matched  filter  or  cross-convolver 
use  of  the  high  capacity  CZT  device  is  shown  in  figure  8.  The  required  column 
multiplexer  may  be  realized  as  a  small  surface  wave  device  with  a  single  input 
tap  and  a  number  of  unweighted  output  taps,  together  with  a  number  uf  high 
speed  analog  switches,  as  shown  in  figure  9.  When  the  device  is  fully  loaded, 
the  switches  are  all  closed  and  the  serially  input  data  column  is  available 
in  parallel  at  the  output  of  the  switches. 


HIGH  CAPACITY  DFT  USING  A  LEXICAL  SCAN 


The  DFT  device  of  figure  7  using  the  linear  congruential  scan  of  equation 
(3)  has  several  Limitations:  the  horizontal  and  vertical  sizes  of  N,;  and 
must  have  no  common  divider,  the  data  must  recirculate  through  the  multiplexer's 
acoustic  surface  wave  device,  and  the  output  of  the  DFT  is  in  the  scrambled 
order  given  by  equation  (4). 

A  lexical  scan  similar  to  one  half  of  a  television  Interlaced  scan  may  be 
used  to  convert  a  one-dimensional  DFT  almost  into  a  two-dimensional  DFT  [7j. 

The  word  "almost"  is  used  becuase  of  the  presence  of  an  additional  phase  factor 
which  must  be  inserted  after  the  first  partial  transform  is  performed.  Although 
this  scan  has  been  previously  proposed  for  parallel  computation  of  the  FFT  [7], 
it  lends  itself  equally  well  to  the  organization  of  parallel  CZT  hardware  and 
largely  circumvents  the  above  mentioned  limitations  of  the  linear  congruential 
scan.  In  particular,  the  lexical  scan  permits  the  complete  elimination  of  the 
large  surface  wave  delay  line  in  the  input  multiplexer. 

In  the  lexical  scan  defined  by  equations  (9)  and  (10),  N  =  niay  be 

any  factorization  of  N  into  the  product  of  two  integers. 


p  =  J1N2  +  J2 


(9) 


q  =  kl  +  k2Nl 


for  *  0,  1 . Nj-1 

J2,k2  *  0,  1,  .  .  .,  N  -1 


(10) 
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The  DFT  of  equation  (1)  may  be  rewritten  in  terms  of  the  lexical  scans 
of  equations  (9)  and  (10)  as  shown  in  equation  (11)  and  simplified  slightly  in 


equation  (12) 


VI V1 


k1+k2N1 


13  Z3fWJ; 

Ji*o  J2-o 


TT  [(h“2  +  V  (ki  +  k2Nl>] 


Lk1+k2N1 


V1  -  ^ 

T,e  e  f3fJivj2  ' 

j2=0  V0 


for  k^  «  0,  1 . N^-l 

k2  “  0|  1)  •  •  •  i  N2”l 


The  structure  of  a  CZT  architecture  using  equation  (12)  is  identical  to 
that  shown  in  figure  7,  but  the  scan  differs  from  that  of  figure  1,  and  the 
required  partial  CZT  differs  from  that  of  figure  4.  The  postmultiplier  of  the 
partial  CZT  for  the  lexical  scan  crn  not  be  the  same  for  all  horizontal  channels 

j2kl 

because  of  the  presence  of  the  additional  factor  of  "*1271  ^  •  In  essence, 

e 

the  postmultiplier  chirps  of  the  partial  CZT  using  the  lexical  scan  must  be  on 
different  (discrete)  carrier  frequencies. 

It  will  be  noted  that  when  a  lexical  scan  is  used,  the  order  in  which  the 
horizontal  and  vertical  transforms  are  performed  is  critical.  For  the  scans 
shown  in  figure  10,  the  vertical  transform  must  be  performed  first. 
Unfortunately,  this  prevents  the  transform  points  from  coming  out  in  natural 
order.  However,  this  should  not  present  any  difficulty  for  performing 
convolution  by  multiplying  in  the  frequency  domain. 


REFERENCES 


[11  Whitehouse,  H.  J.,  J.  M.  Speiser,  and  R.  W.  Means,  High  Speed  Serial 
Access  Linear  Transform  Implementations,  presented  at  the  All  Applications 
Digital  Computer  (AADC)  Symposium,  Orlando,  Florida,  23-25  January  1973, 
reprinted  as  NUC  TN  1026. 

[2 J  Means,  R.  W.,  D.  D.  Buss,  and  H.  J.  Whitehouse,  Real  Time  Discrete 
Fourier  Transforms  Using  Charge  Transfer  Devices,  Proceedings  of  the  CCD 
Applications  Conference  held  at  the  Naval  Electronics  Laboratory  Center,  San 
Diego,  Calif.,  18-20  Sept.  1973,  pp.  95-101. 

[3]  A Is up,  J.  M. ,  R.  W.  Means,  and  H.  J.  Whitehouse,  Real-Time  Discrete 
Fourier  Transforms  Using  Surface  Acoustic  Wave  Devices,  Proceedings  of  the  I EE 
International  Specialist  Seminar  on  Component  Performance  and  Systems 
Applications  of  Surface  Acoustic  Wave  Devices,  held  at  Aviemore,  Scotland, 

24-28  Sept.  1973. 

[4]  Preisendorfer,  R.  W.,  Introduction  to  Fast  Fourier  Transforms,  Visibility 
Laboratory,  University  of  California,  San  Diego,  Spring  1967. 

[5]  Speiser,  J.  M.  and  H.  J.  Whitehouse,  A  Two  Dimensional  Discrete 
Fourier  Transform  Architecture,  NUC  TN  1221,  17  Oct.  1973. 

[6]  Cooley,  James  W. ,  Peter  A.  W.  Lewis,  and  Peter  D.  Welch,  Historical 
Notes  on  the  Fast  Fourier  Transform,  IEEE  Transactions  on  Audio  and  Electro¬ 
acoustics,  Vol.  AU-15,  pp.  76-79,  June  1967. 

[7]  Gold,  Ben  and  Theodore  Bially,  Parallelism  in  Fast  Fourier  Transform 
Hardware,  IEEE  Transactions  on  Audio  and  Electroacoustics,  Vol.  AU-21,  No.  1, 
pp.  5-16,  Feb.  1973. 


D-10 


EXAMPLE 

LET  N,  I 
N>  4 


FiR.  1 


llliislr.ihon  ol  I  lie  iksircd  w.in  lor  N|  '  N  , 


D-ll 


INPUTS 


DEMULTIPLEXER  AND  CHIRP  MULTIPLIER 
(N,  SEPARATE  INPUT  TAPS) 


DISCRETE  CHIRP  ULTER 

(SINGLE  WEIGHTED.  SUMMED  OUTPUT 

OUTPUT  TRANSDUCER  WITH  2N,-1  TAPS) 


a$  =  eW/N,) 


Tup  weights  and  structure  for  acoustic  surface  wave  combined  demultiplexer,  chirp  multiplier, 
and  discrete  chirp  litter. 


IVnodic  convolution  using  2-1)  I'ZT  (output  in  M.r;tmblvd  order) 


SURFACE  WAVE  TAPPED  DELAY  LINE 
CONTAINING  N,  SAMPLES 


CLOSE  EVERY  N 
SAMPLE  TIMES 


Fig.  9  Column  Multiplexer 


APPENDIX  E 


STRUCTURAL  ORGANIZATION 
FOR  REAL  AND  COMPLEX 
CONVOLUTION  BY  IMAGING  CCDs 


Dr.  J.  W.  Bond 

Naval  Undersea  Center 


STRUCTURAL  ORGANIZATION  FOR  REAL  AND  COMPLEX 
CONVOLUTION  BY  IMAGING  CCDs 

By  J.  W.  Bond 
ABSTRACT 

A  mathematical  model  is  developed  for  imaging  CCDs  with  one  or  two  registers. 

It  is  shown  that  an  imaging  CCD  can  be  characterized  by  one  sequence  associated 
to  the  sensors  and  one  sequence  associated  to  each  of  the  registers.  Since  only 
non-negative  amounts  of  charge  :an  be  accumulated  or  transferred,  the  sequences 
associated  with  the  sensors  and  the  registers  must  be  non-negative  and  real.  Meth¬ 
ods  for  calculating  an  arbitrary  real  or  complex  convolution  using  imaging  CCDs 
are  established  by  using  the  time  domain  to  distinguish  the  real  and  imaginary, 
positive  and  negative  terms  needed  C*  calculate  the  desired  convolutions.  Explicit 
architectures  are  described  for  calculating  both  real  and  complex  convolutions 
using  a  CCD  with  either  two  or  only  one  register.  The  recommended  architectures 
for  calculating  complex  convolutions  are  based  on  a  generalization  of  Tiemann  a 
algorithm  for  calculating  a  complex  convolution  by  calculating  the  real  convolu¬ 
tion  of  extended  sequences.  More  generally,  any  number  of  real  and  complex 
convolutions  can  be  calculated  interwoven  in  time  by  a  device  which  can  be  used 
to  directly  calculate  convolutions  of  sequences  of  non-negative  real  numbers. 

INTRODUCTION 

The  chirp  7.  algorithm  reduces  the  calculation  of  a  Fourier  transform  of  a  sequence  of 
numbers  to  a  premultiplication  by  a  complex  chirp,  followed  by  a  complex  convolution  by 
a  complex  chirp,  followed  by  a  post  multiplication  by  a  complex  chirp.  The  CCD  seems  in 
principle  to  offer  a  possibility  of  sensing  and  performing  the  premultipheation  and  convolu¬ 
tion  of  an  image  necessary  to  spatially  Fourier  transform  an  image.  Premultiplication  would 
be  accomplished  by  masks  over  the  sensors,  and  convolutions  would  be  accomplished  by 
gating  the  accumulated  charge  under  the  sensors  to  registers  which  perform  the  required 
delay  and  sum  operations  through  charge  transfer.  This  note  grew  out  of  a  successful  attempt 
to  show  that  CCDs  could  be  used  to  perform  the  Fourier  Transform  apart  from  additional 
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operations  which  do  not  require  memory  or  postnmltiplication.  The  key  result  described 
in  this  note  is  a  scheme  by  which  a  charge-coupled  device  can  be  used  as  a  quite  general 
real  convolver  even  though  only  non-negative  amounts  of  charge  can  be  accumulated  or 
transferred.  Once  this  scheme  has  been  explained  it  is  easy  to  extend  usual  implementations 
of  Fourier  transforms  to  implementations  by  CCDs. 

The  implementations  described  depend  on  an  essential  way  on  the  fact  that  the  image  to 
be  Fourier  transformed,  sampled  by  the  sensors,  is  described  by  a  non-negative  sequence  of 
real  numbers.  The  signs  involved  in  all  calculations  then  are  independent  of  the  actual  data. 
The  image  will  be  assumed  to  be  varying  so  slowly  that  the  image  can  be  assumed  stationary 
during  the  time  required  to  calculate  the  spatial  Fourier  transform  of  the  image  by  the  CCD 
and  its  ancillary  equipment. 

A  useful  mathematical  description  ot  the  CCD  to  derive  the  architecture  for  implement¬ 
ing  complex  convolutions  follows.  Each  CCD  consists  of  line  of  sensor  elements  and  either 
a  right  and  lett  register  each  consisting  ol  data  cells  or  a  single  register  with  data  cells.  Charge 
accumulates  under  the  sensors  and  is  transferred  from  all  sensors  to  some  predetermined  and 
fixed  data  cells  in  the  right  register  or  from  all  sensors  to  some  predetermined  and  fixed 
data  cells  in  the  lett  register.  Alter  charge  has  been  transferred  to  either  register  it  is 
advanced  from  the  data  cell  it  is  in  to  the  one  below  it  before  the  next  charge  is  transferred 
to  either  register.  Figure  I  illustrates  the  basic  one  register  device  fundamental  to  this  note. 
For  ease  of  presentation,  the  number  ol  data  cells  per  register  and  number  of  sensor  elements 
will  be  taken  the  same  so  that  transfers  of  charge  in  the  device  are  always  to  the  data  cells 

adjacent  to  a  sensor  element.  In  figure  I .  g().  gj . gN_  (  denote  the  charge  accumulated 

in  the  sensor  elements  during  a  lixed  accumulation  time  T  and  are  therefore  non-negative 
real  numbers.  These  real  numbers  can  be  viewed  as  a  non-premultiplied  sample  of  the  image 

on  the  sensors.  The  gQ.  g| . g|^_l  can  be  premultiplied  by  numbers  between  0  and  I  by 

using  masks  to  cut  down  the  illumination  hitting  the  sensor  elements.  The  numbers  aQ. 

a  I . dN-l  between  0  and  I  denote  the  proportion  of  illumination  passed  by  the  masks 

on  the  corresponding  sensor  elements.  The  premultiplied  g0.  gj . g^_l  can  all  be  mul¬ 

tiplied  by  the  same  number  between  0  and  I  by  varying  the  time  charge  is  allowed  to  accu¬ 
mulate  in  all  the  sensors.  The  operation  of  the  device  is  described  by  specifying  the  masks 


REGISTER 
DESCRIBED  BY 


SENSOR 

ELEMENT 


SENSOR 


9N-laN-l 

9N-2aN-2 


9lal 


9oao 


CHARGE 

TRANSFER 


DATA 

CELL 


REGISTER 


OUTPUT 

O(t) 


DIRECTION 

REGISTER 

ADVANCED 


PROPORTION 

CHARGE 

TIME 

TRANSFERRED 

r 

t 

O 

0 

ci 

‘l 

• 

• 

• 

CN-1 

tN-l 

CYCLICALLY 

EXTENDED 


N-l 

THEN  0(tN+k)  =  ^ 
j=0 


9jaj<k 


Figure  1.  Imaging  CCD  as  real  convolver. 

and  the  timing  sequence  tor  charge  transfer.  Immediately,  alter  charge  is  transferred  the 
register  is  advanced,  i.e.  the  charge  in  any  data  cell  is  transferred  to  the  data  cell  below  it. 
The  output  of  such  a  device  is  given  by  the  formula  in  figure  I .  which  is  derived  in  Corollary 
Figure  2  illustrates  how  an  imaging  CCD  witli  4  elements  would  be  used  to  convolve  the 
sequence  g0aQ,  g j a  j ,  g2a2-  g3a3  w'tl’  l*,e  sequence  c0,  c  | ,  c-i,  c^  with  the  gj  non-negative 
and  the  aj  and  Cj  between  0  and  1  for  j  =  0,  1,2,  and  3. 


0  (h)  =  9lalco  +  Vocl 
0(t2)  =  g2a2co  +  glalcl  +  90aoc2 


0(13)  =  93a3c0 +g2a2Cl +g1a1c2  +  g0a0c3 

0(t4)  =  g3a3c1  +g2a2c2  +9lalc3  +  90aoco 

0  O5)  =  g3a3c2  +  g2a2c3  +  glalco  +  goaocl 

°(l6)  =  g3a3c3  +  g2a2co+glalcl +9oaoc2 
°(t7)  =  93a3cQ +g2a2Cl +g1a1c2  +  gQa0c3 


NOTE:  AFTER  FIRST  THREE  TRANSFERS  THE  IMAGING  DEVICE 

CYCLICALLY  CONVOLVES  g0aQ,  g j a j ,  g2a2,  g3a3  WITH  c0,  c^, 


Figure  2.  Imaging  C'CT)  as  a  real  convolver  of  non-negative  sequences. 


An  imaging  CCD  can  only  be  used  to  directly  calculate  convolutions  of  two  sequences 
of  non-negative  real  numbers.  In  order  to  use  such  a  device  to  perform  complex  convolutions 
the  following  mathematical  problem  must  be  solved. 

How  can  G  =  g(a+ib)  *  (c+id)  be  calculated  using  only  non-negative  real  numbers  for 
a,  b,  c,  d,  sequences  of  real  numbers  between  -  I  and  I ,  and  g  a  sequence  of  non-negative 
numbers? 

A  straightforward  mathematical  solution  to  this  problem  is  to  write 


|a|-a  ,  .  /|b|+b 


Icl-c  .  /  Id  |+d  Id  |  d 
2  +  \  2 - 2~ 


Note  that  each  of  the  quantities  |x|±x  is  non-negative  so  that  this  identity  calculates  the 
complex  convolution  in  terms  of  16  convolutions  of  non-negative  sequences.  If  this  mathe¬ 
matical  solution  were  implemented  in  the  most  straightforward  maimer  16  CCDs  would  be 
required  each  with  one  register  with  a  total  of  16  sensor  elements  associated  with  each  gj, 

j  =  0 . N-l.  It  could  be  implemented  in  one  CCD  with  I6N  sensor  elements  and  I6N 

data  cells.  Preferred  implementations  are  based  upon  generalizations  of  a  procedure  dis¬ 
covered  by  Tiemann  to  perform  complex  convolutions  with  a  single  real  convolver.  These 
implementations  require  ION  sensor  elements  and  ION  data  cells  when  an  imaging  CCD  is 
used  with  the  same  number  of  data  cells  and  sensor  elements.  If  a  two  register  CCD  were  used 
then  5N  sensor  elements  and  ION  data  cells  are  required.  In  the  process  of  discovering  the 
recommended  implementations  a  general  theory  for  using  CCDs  to  perform  complex  con¬ 
volutions  vas  developed  and  will  be  presented  in  this  note. 

The  approach  is  to  describe  how  a  single  CCD  can  be  used  to  do  calculations  which 
become  increasingly  complex.  First  it  is  shown  how'  to  do  convolutions  of  sequences 

g0aQ, ....  gN_|  ajNj_  |  and  cQ.  Cj . c^_|  when  the  Uj  and  Cj  are  real  numbers  be,  ween 

-I  and  I  with  a  single  CCD  with  two  registers  and  then  with  one  register.  Next,  it  is  shown 

how  to  do  complex  convolutions  of  sequences  g()  (aQ  +  i b Q ) . g^i  ( ii ^ _ j  +  ibj^_| )  and 

cn  +  id() . CN-1  +  '^N-l  ,v'K'n  aj>  bji  cj’  and  dj  are  real  numbers  between  -  I  and  I. 

All  the  implementations  derived  are  based  on  the  decomposition  of  a  real  number  \ 
into  two  non-negative  real  numbers  x+  and  x  defined  by  x+  =  and  x  =  — 
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Other  implementations  are  possible,  in  fact  all  the  formulas  derived  would  hold  lor  |x| 
replaced  by  any  real  number  at  least  as  large  as  |x|.  In  particular,  when  replaced  by  1  the 
implementations  become  implementations  using  “biased  arithmetic”.  Furthermore,  the  for¬ 
mulas  developed  always  involve  subtracting  x  terms  from  corresponding  x+  terms,  i.e.,  they 
make  use  of  the  fundamental  identity  x+  -  x  =  x,  so  that  replacing  x+  by  1  +  x  and  x  by  1 
would  lead  to  still  other  implementations  involving  “biased  arithmetic”. 

Mathematically  a*c  and  b*c  can  be  calculated  by  convolving  a  new  sequence  a(),  b{),  aj, 

b| . aj^_| ,  bj^_|  with  a  new  sequence  cQ,  o,  cj .  o . CN-I>  °-  ^lls  s'mP'e  observation 

is  the  basis  for  reducing  the  numbers  of  CCDs  required  in  implementing  real  and  complex  con¬ 
volutions.  A  simple  notation  for  describing  these  and  similar  sequences  which  can  be  created 

from  given  sequences  follows.  Let  [a,  b,  c . d|  denote  the  sequence  aQ,  bQ,  cQ . du, 

a | ,  b | ,  C| . dj . aN_|,  bN_,,  cN_, - ,  CN_| . dN_,.  where  a  is  the  sequence 

aQ,  a | . aN_j,  b  the  sequence  bc,  b| . bN_( .  c  the  sequence  cQ.  Cj . cN_| . 

d  the  sequence  dQ,  d  j . ^N-l  •  *n  t,lis  notat‘on  l'ie  w'bove  observation  is 

[a,  b  1  *  [c,  01  =  [a*c,  b*c| 
and  Tiemann’s  algorithm  is 

[a.  b,  0|  *  (c,  +d, -e|  =  [#,  a*d  +  b*c,  -  (a*c  -  b * d > ] , 

where  #  is  of  no  interest.  The  recommended  ii>  ’Mementation  of  a  complex  convolution  by 
a  CCD  with  two  registers  is  based  on  convolving 

[a+,  b+,  a  b”,  a+l 

with 

[c+,  d+,  0,  0.  01  by  one  register,  and 

with 

[0,  0,  0,  c" ,  d  1  by  the  other  register. 

The  recommended  implementations  using  one  CCD  with  a  single  register  is  either  to  convolve 
[a+,  a"  ,  b+,  b  .  0,  0,  0.  0,  0,  0| 

with 
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[c+,  0,  d+,  0,  c  ,0,  d  ,0,  c+l  . 


or 


[a+,  0,  b+,  0,  a* ,  0,  b~,0,  a+,  0] 


with 

]c+,  c",d+,  d",0,  0,  0,  0,  0,0]. 

We  note  that  CCD  imagers  could  be  designed  which  would  be  well  suited  tor  performing 
complex  convolutions.  If  each  sensor  element  were  connected  to  the  first  of  two  data  cells 
then  a  one  register  CCD  with  5N  sensor  elements  could  be  used  to  calculate  the  complex  con¬ 
volution  of  two  sequences  of  length  N.  If  more  data  cells  per  sensor  element  were  available 
and  the  particular  data  cell  among  those  associated  to  each  sensor  element  could  vary  from 
sensor  to  sensor  then  a  device  with  2N  sensor  elements  and  10N  data  cells  could  be  used  to 
perform  the  same  complex  convolution. 

USE  OF  A  SINGLE  CCD  TO  PERFORM  CONVOLUTIONS 

In  this  section,  recommended  architectures  are  described  for  performing  a  real  and  com¬ 
plex  convolutions  using  one  and  two  register  CCDs.  The  use  of  CCDs  to  calculate  ga*c 
where  g  is  a  non-negative  sequence  of  real  numbers  and  a  can  c  are  sequences  of  real  numbers 
will  first  be  described  in  this  section.  Since  convolution  is  linear  in  a  and  c,  there  is  no  loss 
of  generality  in  restricting  a  can  c  to  sequences  of  real  numbers  between  -1  and  1 . 

Note  ga*c  =  g(a+-a  )  *  (c+-c  ),  where  a+,  a' ,  c+,  c~  are  the  sequences  obtained  from 
a  and  c  by  defining  a+  =  ( lajl+aj )/2,  aj  =  (|aj|-aj)/2,  cf  =  ( |Cj|+Cj)/2  and  cj  =  ( |cj|-Cj  )/-,. 
Therefore,  if  the  architecture  is  based  on  this  identity,  it  is  necessary  to  calculate  four  con¬ 
volutions  involving  real  sequences  of  non-negative  numbers:  ga+  *  c+,  ga-  *  c+,  ga+  *  c" 
and  ga  *  c~ . 

If  a  CCD  has  two  registers  ga+  *  e+  and  ga  *  c+  can  be  calculated  in  the  right  register 
and  ga+  *  c  and  ga  *  c“  ean  be  calculated  in  the  left  register.  Figure  3  shows  by  example 
how  one  CCD  with  two  registers  can  be  used  to  calculate  ga*c. 

If  a  CCD  has  only  one  register,  ga+  *  c+,  ga"  *  e+,  ga+  *  c~,  and  ga"  *  c~  can  be  calcu¬ 
lated  in  two  different  ways  as  illustrated  in  Figures  4  and  5.  Note  the  implementations 
illustrated  in  Figures  4  and  5  are  obtained  from  one  another  by  interchanging  a  and  e.  The 
figures  in  this  section  present  the  implementations  for  N-3.  It  is  easy  to  derive  the  imple¬ 
mentations  for  general  N  from  those  illustrated  for  N=3  in  the  figures. 
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LEFT  REGISTER  RIGHT  REGISTER 

CHARGE  TRANSFER  CCD  WITH  CHARGE  TRANSFER 

SEQUENCE  MASK  LEVELS  INDICATED  SEQUENCE 
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1  °L<» 
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O  + 

t6 

*7 

co 

• 

• 

1°r(16)  '  °R(t7)  =  goaoco  +  glalc2  +  g2a2cl] 

-  [0L(t7)-0L(t8)  =  90a0c0+  g1alc2*q2a2c\] 

=  goaoco  +glalc2  +g2a2cl 
[Or (t8)  -  0R(tg)  =  90aocl +9la1co  +  g2a2C2] 

-  (0L(t9)  -  oR(t10)  =  g0a0cI  +g1a1c0  +  g2a2c2] 

=  goaocl  +  glalco  +  g2a2c2 

(Or^K))  OR^ll)  =  goaoc2  +  glalcl  +  g2a2co) 

-  (0L(tn)  -  0R(t12)  =  g0a0c2  +  g1a1c1 +g2a2c0] 

=  goaoc2  +glalcl  +  g2a2co 

Figure  3.  Real  convolution  of  length  3  sequences  with  a  2  register  CCD. 


CCD  WITH 
MASK  LEVELS 
INDICATED 


CHARGE  TRANSFER 
SEQUENCE 


0(tj2)  “  OU13)  "  0(^14)  +  0^15)  =  90aoco  +  Q\*\c2  +  92a2cl 

0(ti6)  -  0(ti7)  -  O(tig)  +  0(t19)  =  90aocl  +  9iai-0  +  92a2c2 

0(t2o)  -  0(t2i)  -  0(^22)  +  0(^2^  =  90aoc2  +  9lalcl  +  S2a2co 


Figure  4.  Rea!  convolution  of  length  3  sequences  with  a  one  register  CCD,  i. 


CCD  WITH  CHARGE  TRANSFER 

MASK  LEVELS  INDICATED  SEQUENCE 
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90ao 
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0 

ho 

|o(t)  0 
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co 
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hi 

t12 

t13 

• 

• 

0(t12)-0{t13)-0{t14) +  0(t15)  =  g0a0c0  +  g1a1c2  +  g2a2cl 

°(t i6>  —  0 1 7)  “  °(t18^  +  °{t19>  =  9oaocl +  9lalco +  92a2a2 
O(t20)  -  0(t21)  -  0{t22)  +  0(t23)  =  g0a0c2  +  gjajCj  S2a2c0 

Figure  5.  Real  convolution  of  length  3  sequences  with  a  one  register  CCD.  II 
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The  implementations  of  the  complex  convolution  g(a+ib)*(c+id)  which  are  described 
next  utilize  Tiemann’s  algorithm.  Note  that  16  sums  would  be  required  if  the  basis  for  the 
calculation  weie  the  identity  g((a+-a_)  +  i(b+-b-))  *  ((c+-c“)  +  i(d+-d-))  =  g(a+ib)  *  (c+id). 
Tiemann’s  algorithm  can  be  written  [a,  b,  0]  *  [c,  d,  -c]  =  [#,  ad+bc,  -(ac-bd)] ,  where  # 
denotes  a  term  which  is  not  useful.  If  implemented  by  a  CCD  in  terms  of  a+-a“,  b+-b“, 
c+-c“,  d+-d“,  1  2  sums  would  be  required.  The  implementations  for  the  complex  convolu¬ 
tion  described  next  are  based  on  a  modification  of  Tiemann’s  identity.  The  two  register 
single  CCD  implementation  is  based  on  the  convolution  of  [a+,  b+,  a-,  b“,  a+]  with  [c+,  d+, 
0,  0,  0]  in  one  register  and  [c“,  d”,  0,  0,  0]  in  the  other.  The  one  register  single  CCD  imple¬ 
mentations  are  based  on  either  the  convolution  of  [a+,  a-,  b+,  b~,  0,  0,  0,  0,  0,  0]  with 
[c+,  0,  d+,  0,  c",  0,  d“,  0,  c+,  0]  or,  by  intercl.  nging  a  and  b  with  c  and  d,  the  convolution 
of  [a+,  0,  b+,  0,  a-,  0,  b~,  0,  a+,  0]  with  [c+,  c“,  d+,  d“,  0,  0, 0, 0,  0,  0] .  Note  these  imple¬ 
mentations  only  require  the  calculation  of  10  sums,  two  of  which  are  not  combined  with 
the  other  input  used  to  form  the  desired  real  and  imaginary  parts  of  the  output  Figure  6 
illustrates  the  2  register  single  CCD  implementation  and  figures  7  and  8  the  two  different 
one  register  single  CCD  implementations  of  the  complex  convolution. 

The  implementations  illustrated  in  figures  4  and  7  have  the  different  masks  of  the 
same  signal  sample  adjacent  while  these  illustrated  in  figures  5  and  8  have  them  separated  by 
sensor  elements  masked  to  zero.  Therefore  for  CCDs  with  each  sensor  element  correspond¬ 
ing  to  a  single  data  cells  the  implementations  described  by  figures  4  and  7  seem  best  because 
adjacent  sensor  elements  could  be  located  within  one  circle  of  resolution  of  the  imager. 

Note  that  in  the  CCD  imagers  illustrated  the  basic  timing  sequence  for  charge  transfers 
can  be  longer  than  the  number  of  sensor  elements  because  zeros  at  either  the  beginning  or 
end  of  the  sequence  associated  to  the  sensors  need  not  have  corresponding  sensor  elements. 
In  particular,  note  that  the  implementation  of  the  complex  convolution  illustrated  by  figure 
7  for  the  case  N  =  3  would  in  general  require  ION-6  sensor  elements  while  that  illustrated  by 
figure  8  would  require  ION-1 

The  difference  in  the  total  number  of  sensor  elements  required  for  the  implementa¬ 
tions  illustrated  in  figures  7  and  8  is  greatly  exaggerated  by  taking  N  =  3  and  this  difference 
is  not  considered  to  be  as  important  as  possible  design  modifications  which  can  be  used  to 
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Figure  8.  Complex  convolution  by  a  one  register  CCD  of  sequences  oi  length  3, 11. 


reduce  the  number  of  sensor  elements  proportional  to  N  for  some  implementations  and  not 
others. 

Suppose  that  a  single  register  CCD  were  designed  so  that  each  sensor  element  trans¬ 
ferred  its  charge  to  the  first  of  two  data  cells.  Then  figure  5  illustrates  an  implementation 
of  the  real  convolution  of  sequences  of  length  N  which  could  be  accomplished  with  a  CCD 
of  this  type  with  2N  sensor  elements  and  4N  data  cells  and  figure  8  illustrates  an  implemen¬ 
tation  of  the  complex  convolution  of  sequences  of  length  N  which  could  be  accomplished 
with  5N  sensor  elements  and  ION  data  cells. 

Suppose  the  number  of  data  cells  in  the  CCD  could  differ  from  the  number  of  sensor 
elements  and  furthermore  that  the  particular  data  cell  to  which  charge  is  transferred  from  a 
sensor  element  could  depend  on  the  location  of  the  sensor  element.  A  two  register  CCD 
designed  in  this  manner  could  convolve  two  real  sequences  of  length  N  using  N  sensor  ele¬ 
ments  and  4N  data  cells.  This  is  because  for  any  real  number  x  either  x+  or  x“  will  be  zero 
and  therefore  the  sign  can  be  taken  into  account  by  the  data  cell  to  which  charge  is  trans¬ 
ferred.  Figure  9  illustrates  a  scheme  which  would  allow  the  implementation  illustrated  in 
figure  3  to  be  accomplished  with  only  2N  sensor  elements.  Figure  10  provides  a  scheme  for 
accomplishir  g  the  complex  convolution  illustrated  in  figure  7  using  2N  sensor  elements  and 
ION  data  cells  with  such  a  one  register  CCD. 

The  2  register  implementation  of  complex  convolution  by  such  CCDs  is  not  so  clear 
cut,  because  the  Tiemann  algorithm  is  not  symmetrical  in  a,  b  and  c.  d.  The  implementa¬ 
tions  with  a  two  register  CCD  illustrated  in  figure  6  has  about  5/2  N  sensors  not  masked  to 
zero  and  uses  ION  data  cells.  The  implementation  of  the  Tiemann  algorithm  based  on 
Corollary  9  for  a  two  register  CCD  has  2N  sensors  not  masked  to  zero  and  would  require 
12N  data  cells.  In  both  cases  the  scheme  in  figure  10  could  be  used  to  do  the  implementa¬ 
tions  using  a  minimal  number  of  sensor  elements  for  the  given  number  of  data  cells. 


E-l  6 


CORRESPONDING 

DATA  CELLS  IN  SENSOR 

LEFT  REGISTER  ELEMENT 


CORRESPONDING 
DATA  CELLS  IN 
RIGHT  REGISTER 


«i  'V 

- 

— 

Figure  9.  A  scheme  for  performing  a  real  convolution  with  a  two  register  CCD 
using  a  minimal  number  of  sensor  elements. 


E-17 


aj  >  0  bj  >  0 


3j>0  bj  <  0 


□ 

I'i1 

Ibjl 

-□ 

□ 

j  <  0  bj  >  0 

The  one  register  connections  would  be  to  the  first  4  data  cells  as  illustrated,  but 
to  only  one  side,  and  the  number  of  data  cells  per  pair  of  sensor  element  would  be 
10. 


Figure  10.  Scheme  for  CCDs  to  implement  complex  convolution  with 
minimal  numbers  of  sensor  elements. 


1.  THE  MATHEMATICAL  THEORY  OF  CALCULATION  WITH 
ONE  AND  TWO  REGISTER  CCDs. 

In  this  section,  we  develop  a  simple  way  to  describe  the  output  of  a  single  and  double 
register  CCD  as  a  function  of  time.  This  will  enable  us  to  derive  the  results  discussed  in  the 
introduction  and  previous  sections  of  this  note  in  a  simple  manner. 

Suppose  a  CCD  consists  of  N  sensor  elements  and  a  single  register  consisting  of  N 
data  cells.  Suppose  further  that  gQ . g^_l  denotes  the  charge  accumulated  by  the  N  sen¬ 

sor  elements  in  a  fixed  accumulation  time  T.  Suppose  that  by  masking  a  sensor  element, 
the  amount  of  charge  accumulated  by  that  sensor  element  can  be  reduced.  Thus  masking 
allows  us  to  multiply  the  non-negative  gj  by  aj,  real  numbers  between  0  and  1 .  The  charge 
accumulated  by  the  j-th  sensor  element  is  written  gjaj;  the  gj  is  explicitly  included  in  some 
of  the  propositions  and  corollaries  because  in  applications  gj  will  depend  on  time  and  aj  will 
not.  Thus  the  architecture  to  implement  a  calculation  will  not  depend  on  the  gj  but  will 
depend  on  the  aj. 

To  describe  the  output  of  the  CCD  let  cQ,  cj . cm, .  .  .  denote  the  proportion 

of  charge  which  could  be  accumulated  during  fixed  accumulation  time  T  which  is  trans¬ 
ferred  from  each  sensor  element  to  its  corresponding  data  cell  at  times  tQ,  tj . tm.  .  . . 

respectively.  Without  loss  of  generality  all  the  Cj  can  be  taken  to  be  real  numbers  between 
0  and  1 .  The  output  of  the  device  at  one  of  these  times  tm  will  be  denoted  O(m)  and  is  the 
amount  of  charge  accumulated  in  the  data  cell  corresponding  to  the  sensor  element  indexed 
by  0.  Immediately  after  a  charge  transfer  the  register  is  advanced  and  the  charge  accumulat¬ 
ed  in  a  data  cell  is  transferred  to  the  next  lowest  indexed  data  cell.  The  next  proposition 
describes  O(ni).  See  figures  1  and  2. 

PROPOSITION  1 

The  output  of  a  CCD  with  N  sensor  elements  and  N  corresponding  data  cells  in  its 
register  is  given  by: 

O(ni)  -  (cm,  cm_i, .  .  .  .  cm_(N-|))  ’  ao’ ai . &N-I  ;1N-P 


f 


where 


1  g0a0,  g,aj , . . .  ,  g]s|-l  aN-l  denotes  the  char&e  accumulated  in  the  accumulation 
period  T 

-i  r,  c  denotes  the  proportion  of  charge  transferred  from  all  sensor 

—  >  w2’  •  •  ■  ’  nv  ‘  ‘ ' 

elements  to  corresponding  data  cells  at  times  tQ) . . .  ,  tm, . .  . 

3.  O(m)  denotes  the  output  of  the  CCD  at  time  tm 

4.  cj  =  0  for  j  <  0  by  convention 

Proof 

The  proof  will  be  by  induction  on  m  and  N.  If  N=1  the  output  formula  is  clearly  true 
for  all  m.  Suppose  that  the  output  formula  has  been  established  for  CCD  with  less  than  N 
sensor  elements.  Suppose  we  have  a  CCD  with  N  sensors.  Then  0(0)  -  (c0,  0, ....  0). 

(g0  a0, .  . .  ,  gN_j  aN_j)  =  g0  aQ  cQ  so  that  the  output  formula  is  true  at  time  tQ.  Suppose 

the  outputs  for  times  tm  <  tm+j  are  given  by 

O(m)  =  (cm,  cm_  j i  ))•  (g0  aQ, .  .  .  ,  gpg_  j  a^.  j ) . 

Now  0(m+1 )  =  cm+1  gQ  aQ  +  Cm,  wh  *re  Cm  is  the  charge  in  the  data  cell  correspond¬ 
ing  to  g0  aQ  after  the  charge  has  been  transferred  at  time  tm  and  has  been  shifted  down.  Cm 

is  therefore  the  output  at  time  tm  of  the  CCD  with  N-l  sensors  associated  with  gj  a, . 

gN_i  aN_, ,  and  the  same  charge  transfers.  Therefore  by  induction  hypothesis  on  N  we  have 

cm  =  tcm'cm-l*-*-*cm-(N-2))*(8l  al  >  -  -  •  ^N-l  aN-P 

=  ^°>  c(m+l)-l’  c(m+l)-2»  •  •  •  ’  C(m+1)-(N-1))  ’  {go  ao>  ■  ■  •  •  gN-l  “N-P 

Therefore 

°(m+l)  =  (cm+l’c(m+l)-l’---’c(m+l)-(N-l))'(goao’'-->gN-l  aN-P 
which  completes  the  proof. 

Note  that  the  sensors  have  been  indexed  from  the  bottom  up  while  charge  is  trans¬ 
ferred  from  the  top  down.  This  is  necessary  so  that  the  imaging  CCD  can  be  viewed  natural¬ 
ly  as  a  real  convolver. 
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Corollary  1 

The  output  of  a  CCD,  when  cm  =  ck  if  m^k  mod  N  with  o  <  k  <  N-l ,  is  given  by 
N-l 

O  ((N)  +  k)  =  ^  gjaj  c'k.j,  for  o  <  k  <  N-l  . 
j=0 


Proof 


0(N+k)  -  (cN+k,cN+k_,,.  ..  ,ck)  •(g0a0,g1a1,...,gN_|  aN_,). 


i  ^  i 

Gk  =  ^  gjaj  ck_j,  k=0, .  .  . ,  N-l 
j=0 

then  G  -  (ga)  *  c.  Thus  Corollary  1  shows  that  the  CCD  can  be  viewed  as  a  device  which 
convolves  the  sequence  (gQa0 - -  gN_,  aN_, )  with  the  sequence  (c0, _ cN_,).  Fur¬ 

thermore  the  CCD  is  completely  described  by  specifying  the  two  sequences. 

For  sequences  a,  b,  c, . .  .  ,  d  all  of  length  N  let  (a,  b,  c, .  .  .  ,  d]  denote  the  sequence 

ao<  bo’  co>  •  •  •  ’  do>  al  ’  b  1 .  C1  ’  •  •  •  .  d  1 . aN-l>  bN-l '  CN-1  >  •  •  •  ’  dN-l  • 

PROPOSITION  2 

Let  a,  b,  c  be  sequences  of  length  N.  Then 
[a,  b]  *  [c,  o]  =  [a  *c,  b  *  c)  . 


Proof 


Let  d  [a,  b]  so  that  d2k  -  ak  and  ^2k+l  =  bk  and  ^  =  °1  so  tbat  ^2k  =  ck  and 


f2k+l  =0-  The" 


^1N-I  IN- 1  N-l 

°2k  =  £  dj'  f2k-j'  =  £  d2j  f2k-2j  =  £ 


ak  ck-j 


(d*f>2k+l  =  J  dj'*(2k+lH'  =  J  d2j+l  f(2k+l )-(2j+l)  J  bkck-j 
j'=0  J=o  j=0 


If  x  is  any  real  number  let 

x+  =  (|x|  +  x)/2  and  x"  =  (|x|  -  x)/2. 


Corollary  2 

A  CCD  with  one  register  can  be  used  to  convolve  ga  with  c,  with  a  any  sequence  of  N 
real  numbers  between  -1  and  1  and  c  a  sequence  of  N  non-negative  numbers  by  summing 
pairs  of  successive  outputs. 

Proof 

Let  the  CCD  be  the  one  associated  with  the  sequences  [ga+,  ga"  1  and  [c,  o]  and  note 
that  (ga+  *  e)  -  (ga-  *  c)  =  g(a+  -  a-)  *  c  =  ga*c. 

Corollary  3 

A  single  CCD  with  2  registers  can  be  used  to  convolve  ga  with  c,  where  a  and  c  are 
sequences  of  real  numbers  between  -1  and  1  by  summing  4  successive  outputs  two  from 
one  register  and  two  from  the  other. 

Proof 

Let  the  sensor  elements  be  associated  with  the  sequence  [ga+,  ga-  ]  and  one  register 
with  the  sequence  [c+,  o]  and  the  other  register  with  the  sequence  [o,  c"  1 .  Then 
(ga+  -  ga- )  *  (c+  -  c“ )  =  ga  *  c. 

Note  by  substituting  for  a,  b,  c  in  Proposition  2  the  sequences  [c,  o] ,  [d,  o] ,  and 
[a,  b],  respectively,  it  follows  that  [(a,  b] ,  (o,  o]]  *  [[c,  o] ,  [d,  o]]  =  [[a,b]  *[c,o), 

(a,  b]  *  [d,  ol]  =  [a  *  c,  a  *  d,  b  *  c,  b  *  dl .  It  is  easy  to  generalize  this  observation  to 
obtain  the  next  theorem. 

Theorem  1 

If  a,  b,  c . d  are  m  sequences  of  length  N  and  e,  f . g  are  n  sequences  of 

length  N  then 

[[a,  b,  c, .  .  . ,  d  1 ,  [o . o], . . . ,  [o . ol]  * 

[[e,  o . o ] ,  [ f ,  o . o] . [g . ol] 

=  la  *  e,  a  *  f . a*g,  b*e,b*f . b  *  g . d  *  e,  d  *  f, .  .  .  ,  d  *  g] 


Proof 

First  we  establish  the  theorem  for  the  case  when  m  is  arbitrary  and  n=l.  Let 

g=  [a,b, . . .  ,d]  and  h  =  [f, . ..  ,o).  Then  gmk  =  ak>gmk+l  =  bk>  ■  ■  •  >  gmk+(m-l ) =  dk 
and  hmk  =  t'k  with  hmk+j  =  0  for  j  =  1, . . . ,  m  -  1.  Therefore 

mN-1 

(8  *  11  W+k'  =  ^  g(mk+k')-j'  hj' 

j“o 

N-l 

-  £  *(mk+k')-mikmj.fOThj'  =  0unlcssj'  =  mj 
j=0 


N-l 

=  ^  gm(k-j)+k'  hmj  • 
j=0 


Therefore 


=  N-l 

N-l 

r“* 

(g  *  hW  =  ^ 

ak-j-fj’  *  h^mk+l  “ 

j=o 

j=o 

n-l 

-I 

dk-j  fj 

j=o 

bk_j  fj»  •  •  •  » (s  *  k)mk+m-l 


completing  the  first  part  of  the  proof. 

Now  the  proof  is  completed  by  induction  on  n+m.  Suppose  the  proposition  has  been 
proved  for  n+m  <  k.  Suppose  n+m  =  k  and  in  addition  that  n  3*  2  for  the  induction  step 
has  been  proved  already  when  n  =  1. 


Let 

A  =  [a,  b,  c, . . .  ,  d] ,  O  =  [o . o),E  =  [e,o,...,o], 

F  =  [f,o,.  .  .  ,o] . G  =  [g . oj. 

Then 


(a,  b,  c, .  .  .  ,  d  1 ,  (o . o), - lo - ,ol  *  [[e.  o . ol,(f.o . o], 

. (g,  o, . . .  ol]  =  (A,0 . 01  ME,F . G1  =  [A  *  E,  A  *  F, ....  A  *  G] 

=  [a  *  e,  a  *  f . a  ♦  g,  b  *  e,  b  *  f . b  *  g . d  *  e.  d  *  f . d  *  g) 


This  completes  the  proof  of  the  theorem. 
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Corollary  4 

A  single  CCD  can  be  used  to  calculate  any  number  of  real  convolutions  interwoven 
in  time  provided  the  sequences  involved  are  non-negative  sequences  of  real  numbers 
bounded  by  1 . 

Corollary  5 

A  single  CCD  with  one  register  can  be  used  to  calculate  a  convolution  of  a  sequence 
ga  and  c  with  a  and  c  both  sequences  of  real  numbers  between  -1  and  1  by  summing  4  suc¬ 
cessive  outputs  of  the  register.  In  particular,  let  the  sensors  be  associated  with  the  sequence 
[a+,  a-,  o,  o]  and  the  register  associated  with  the  sequence  [c+,  o,  c~,  ol . 

Proof 

[a+,  a",  o,o]  *  [c+,o,  c",o|  =  (a+  *  c+,  a'  *  c+,  a+  *  c~ ,  a”  *  c'  ] 
and 

(a+  *  c+)  -  (a+  *  c“)  -  (a-  *  c+)  +  (a-  *  c~)  =  a  *  c. 

Corollary  6 

One  CCD  with  two  registers  can  be  used  to  calculate  a  complex  convolution 
(ga  +  ib)  *  (c+id)  in  a  conventional  manner  with  a,  b,  c,  and  d  real  numbers  between  -1 
and  1.  In  particular,  let  the  sensors  associated  with  one  CCD  be  [a+,  a",  b+,  b~,  o,  o,  o,  o] 
with  one  register  [c+,  o,  o,  o,  d+,  o,  o,  o]  and  the  other  register  [o,  o,  c”.  o,  o,  o,  d~ ,  o] . 

Proof 

The  right  register  calculates 

[a+,a~,b+,  b‘,o,  o,  o,  ol  *  [c+,  o,  o,  o,  d+,  o,  o,  o) 

=  (a+  *  c+,  a-  *  c+,  b+  *  c+,  b"  *  c+,  d+  *  d+,  a  *  d+,  b+  *  d+.  b  *  d+] 

so  that 

1 )  (a+  *  c+  -  a'  *  c+)  -  (b+  *  d+  -  b-  *  d+)  =  a  *  c+  -  b  *  d+ 
and 

2)  (b+  *  c+  +  b~  *  c+)  +  (a+  *  d+  -  a'  *  d+)  =  b  *  c+  +  a  *  d+. 
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The  left  register  yields  the  same  convolutions  as  the  right  register  with  e+  leplaeed  by 
e  and  d+  replaced  by  d~  so  that 

(a  *  c+  -  b  *  d+)  -  (a  *  e"  -  b  *  d' )  =  a  *  c  -  b  *  d 
and 

(b  *  c+  +  a  *  d+)  -  (b  *  c“  +  a  *  d" )  =  b  *  e  +  a  *  d. 

Corollary  7 

One  CCD  with  one  register  can  be  used  to  calculate  a  complex  convolution  g  (a+ib)  * 
(e+id)  in  a  conventional  manner  with  a,  b,  c,  and  d  real  numbers  between  -  I  and  1.  In 
particular,  let  the  sensors  be  associated  with  the  sequence 
[a+,  a",  b+,  b  ,  o,  o,  o,  o,  o,  o,  o,  o,  o,  o,  o,  ol 
and  the  register  with 

[e+,  o,  o,  o,  e" ,  o,  o,  o,  d+,  o,  o,  o,  d~ ,  o,  o,  o  ] . 

Proof 

The  convolution  ot  the  two  sequences  clearly  gives  the  same  sequences  a+,  a-,  b+,  b- 
convolved  with  c+,  c_,  d+,  d~  as  used  above  to  construct  a  *  c  -  b  *  d  and  b  *  c  +  a  *  d. 

Corollaries  6  and  7  do  not  yield  the  best  implementation  of  a  complex  convolution 
using  a  CCD,  for  a  better  one  can  be  derived  based  on  an  algorithm  found  by  Tiemann. 
Tiemann’s  procedure  is  an  immediate  corollary  to  the  next  proposition. 

PROPOSITION  3 

Let  a,  b,  c,  d,  and  e  be  sequences  of  N  real  numbers.  Then 
[a,  b,  o)  *  [c,  d,e|  =  [#,  a  *  d  +  h  *  c,  a  *  e  +  b  *  d) 
where  the  expression  for  #  in  terms  of  a,  b,  c.  d,  and  e  is  not  of  interest. 

Proof 

Let  f  =  [a,  b,  c]  and  g  =  [c,d,  el  so  that  f3k  =  ak,  f3k+)  =  bk,  l'3k+2  =  0,g3k  =  ck, 
g3k+l  =  dk’andg3k+2  =  ek-  Tllen 


1 
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N-l  N-l 

(f*g)3k+l  =  2  f3jg(3k+i)-3j  +  2  f3j  g(3k+l )-(3j+l ) 

j=0  j=0 

N-l  N-l 

=  I  ak  dk-j  +  2  bk  ck_j 

j=0  j=0 

N-l  N-l 

(f*8)3k+2  =  2  f3j  g(3k+2)-3j  +  2  f3j+l  g(3k+2)-(3j+l ) 
j=0  j=0 

N-l  N-l 

=  2  ai  ek-j +  2  bj  dk-j 

j=0  j=0 

Corollary  8  (Tiemann) 

Let  a,  b,  c  and  d  be  sequences  of  N  real  numbers.  Then 
[a,  b,  o]  *  (c,  d,  -cl  =  [#,  a*d+b*c. -(a*c-b*  d)] 
where  the  expression  for  #  in  terms  of  a,  b,  c  and  d  is  of  little  interest. 

Proof 

Let  e  =  -c  in  proposition  3. 

Corollary  9 

One  CCD  with  two  registers  can  be  used  to  implement  the  Tiemann  algorithm  for 
complex  convolutions  In  particular,  let  the  sensors  be  associated  with  the  sequence 
[a+,  a-,  b+,  b”,  o,  o]  and  one  register  with  the  sequence  [c+,  o,  d+,  o,  c~,  o]  and  the 
other  register  with  the  sequence  (o,  c” ,  o,  d“ ,  o,  c+) . 

Proof 

Let  A  =  [a+,  b+,  o] ,  B  =  [a-,  b",  o),C  =  (c+,  d+,  c"l  and  O  =  (o,  o,  o) 
in  proposition  3.  Then 

[a+,  a",  b+,  b”,  o,  o]  *  [c+,  o,  d+,  o,  c' ,  o) 

=  (A,  B|  *  [C,0]  =  [A*C,  B*C1 
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Note 

A*C  =  |#,  a+  *  d+  +  b+  *  c+.  (a+  *  e"  +  b+  *  d+)l 
and 

B*C  =  I#,  a-  *  d+  +  b_  *  c+,  (a-  *  c"  +  b“  *  d+)  by  Proposition  3. 

Then 

a+  *  d+  +  b+  *  c+  -  (a~  *  d+  +  b"  *  c+)  =  a  *  d+  +  b  *  c+ 
and 

(a+  *  c~  +  b+  *  d+)  -  (a“  *  c“  +  b-  *  d+)  =  a  *  c"  +  b  *  d+  . 

The  left  register  outputs  can  be  used  to  calculate  the  same  expressions,  but  with  e+ 
replaced  by  c~,  d+  by  d~,  and  c~  by  c+.  There! ore 

(a  *  d+  +  b  *  c+)  -  (a  *  d~  +  b  *  c“)  =  a  *  d  +  b  *  c 

and 

(a  *  c+  +  b  *  d~)  -  (a  *  c~  +  b  *  d+)  =  a  *  e  -  b  *  d  . 

Corollary  10 

One  CCD  with  one  register  can  be  used  to  implement  the  Tiemann  algorithm  for  com¬ 
plex  convolutions.  In  particular,  let  the  sensors  be  associated  with  the  sequence  |a+,  a  ,  b+. 
b~,  o,  o,  o,  o,  o,  o,  o.  o]  and  the  register  with  the  sequence  |c+,  o.  d+.  o.  c~.  o.  e“.  o.  d",  o. 
c+,  o] . 

Proof 

[|a+,  b+,  o,  o,  o.  ol .  ia~.  b".  o.  o,  o,  ol]  *  [| .  d+.  c".  c".  d".  c+l  *  |o.  o.  o,  o.  o.  o.l] 

yields  all  the  outputs  used  to  construct  a  *c  -  b  *d  and  a  *  d  +  b  *  c  in  the  Corollary  9. 

Note  that  the  register  sequence  described  in  Corollary  10  repeats  c  in  the  middle. 

This  suggests  that  length  ION  sequences  might  be  used.  Towards  this  end  the  next  proposi¬ 
tion  is  established. 
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PROPOSITION  4 


Let  a,  b,  c,  d,  e,  f,  and  g  denote  sequences  of  length  N.  Then 

[a,  b,  o,  o,  o]  *  [c,  d,  e,  f,  g]  =  [#,  a*d  +  b*c,  a*e  +  b*d,  a*f  +  b*e,  a*g  +  b*f) 
where  the  expression  for#  is  not  of  interest. 


Proof 

Let  p  =  [a,  b,  o,  o,  o]  and  q  =  [c,  d,  e,  f,  g] .  Then  p5k  =  ak,  p5k+j  =  bk,  p5k+j  =  0, 
j  =  2.  3,  4  and  q5k  =  ck,  q5k+1  =  dR,  q5k+2  =  ek,  q5k+3  =  fk,  and  q5k+4  =  gk.  Then 


N-l  N- 1 

(P*q>5k+h  =  ^  P5j  q(5k+h)-5j +  J  p5j+l  q(5k+h)  -  (5j+l) 
j=0  j=0 


N-l  N-l 

=  S  p5jq5(k-j)+h+  £  P5j+I  q5(k-j)  +  (h-l)  -  h  =  0,l,2.3. 
i=n  i=n  and  4. 


Thus  (p  *  q)5k+h  gives  the  desired  terms  for  h  =  1,2,  3,  and  4. 


Corollary  1 1 

One  CCD  with  two  registers  can  be  used  to  implement  a  modifed  Tiemann  algorithm 
for  complex  convolutions.  In  particular,  let  the  sensors  be  associated  with  [a+,  b+,  a“,  b“, 
and  one  register  with  [c+,  d+,  0,  0,  01  and  the  other  with  [c“,  d“,  0,  0,  0] . 


Proof 

[a+,  b+,  a-,  b",  a+]  *  [c+,  d+,  0,  0,  0]  =  [#,  b+  *  e+  +  a+  *  d+,  a'  *  c+  +  b+  *  d+, 
b~  *  c+  +  a”  *  d+,  a+  *  c+  +  b-  *  d+] .  Then  (b+*c+  +  a+*d+)  -  (b”*c+  +  a~*d+)  =  b*c+ 

+  a*d+  (a+*c+  +  b“*d+)  -  (a~*c+  +  b+*d+)  =  a*c+  -  b*d+.  The  output  of  other  register 
can  be  used  to  obtain  the  same  expressions  with  c+  replaced  by  c-  and  d+  replaced  by  d~. 

Corollary  12 

One  CCD  with  one  register  can  be  used  to  implement  a  modified  Tiemann  algorithm 
for  complex  convolutions.  In  particular,  let  the  sensors  we  associated  with  [a+,  a",  b+,  b“, 
o,  o,  o,  o,  o,  o]  and  the  register  with  the  sequence  [c+.  o,  d+,  o,  c~,  o,  d",  o,  c+,  ol . 


Proof 


Let  A  =  (a+.  b+,  o,  o,  o],  B  =  [a“,  b“  o,  o,  o] ,  C  =  (c+,  d+,  c',  d”,  e+]  and 
O  =  [  o,  o,  o,  o,  o] .  Then 

[A,  B]  *  [C,0]  =  [A*C,  B*C] . 

As  in  the  proof  of  Corollary  9,  the  expressions  a*d+  +  b*c+,  a*c~  +  b*d+,  a*d“  +  b*e“,  and 
a*c+  +  b*d+  are  obtained  by  taking  the  difference  of  successive  terms.  Then 

(a*d+ +  b*c+)  -  (a*d“+ b*c~)  =  a*d  +  b*c  and  (a*c+ +  b+d“)  -  (a*e“  +  b*d+) 

=  a*c  -  b*d  . 

Because  of  the  importance  of  Corollary  1 2  it  is  desirable  to  have  a  direct  proof  of  it 
independent  of  most  of  the  theoretical  development  presented  in  this  section. 


Direct  Proof  of  Corollary  12 
Let 

e  =  [a+,  a",  b+,  b~  o,  o,  o,  o,  o,  o] 
f  =  [c+,o,d+,  o,  c“,o,d“,  o,  c+,  o] 

Then 

e  1  Ok  =  ak+’e10k+l  =  ak~  e10k+2  =  bk+’e10k+3  =  bk~  e10k+j  =  0  for 
j  *  4, ....  9 
and 

f  1  Ok  =  ck+’  f10k+2  =  dk+>f10k+4  =  ck~’f10k+6  =  dk~’f10k+8  =  ck+  and 
f  1  Ok+j  =  0  when  J  =  L3,5,  7,  and  9. 

The  proof  will  follow  by  calculating  the  8  useful  terms  of  e  *  f. 

First 


Nr.1  N-l 

(e*f)10k+h  2  e10j  f(IOk+h)  -10j  +  ^  e  lOj+1  f(10k+h)  -  ( lOjf  1 ) 

j=0  j=0 

N^l  N-l 

+  2  e10j+2f(40k+h)-(10j+2)+  2  e  1  Oj+3  *(  1  Ok+li)  -  ( 1  Oj+3) 
j=0  j=0 
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because  e|0k  +  j*  =  0  for  k'  =  4, . . . ,  9.  Next 

N-l  N-l 

(e*0 1  ok  +  2h  =  ^  e10jf10(k-j)  +  2h  e10j+2  f10(k-j)  +  2(h-l) 
j=0  j=0 

and 


N-l  N-l 

(e*^10k  +  (2h+l)  ^  e10j+l  f10(k-j)  +  2h  e10j+3  f10(k-j)  +  2(h-l ) 

j=0  j=0 

because  f]  ok  +  k'  =  ®  when  k'  is  odd. 

Therefore 


N-l 

gk  *  (e*f)10k  +  2h"(e*^10k  +  (2h+l)  =  J  (e  1  Oj  ~e10j+l)f10(k-j)  +  2h 

j=0 

N-l 

+  S  (e10j+2-e10j+3)f10(k-j)  +  2(h-l)- 
j=0 

Then 


N-l 

N-l 

g,  =  j 

ai  dk-i +  £ 

j=o 

j-0 

N-l 

N-l 
■  » 

g2  =  j 

aj  ck-j  +  z 

bj  dk-j 

j=o 

j=o 

N-l 

N-l 

g3  =  £ 

aj  dk-j +  z 

bj  ck-j 

j=o 

j=0 

N-l 

N-l 

a 

J 

aicN+  S 

bj  dk-j 

j=o 

j=0 

so  that  finally 
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ABSTRACT 


This  note  describes  additions  and  modifications  to  a  signal  processing 
interpreter  which  was  developed  as  part  of  the  ARPA  sponsored  program  for 
image  transmission  via  spread  spectrum  links.  The  signal  processing  inter¬ 
preter,  SPIN3, is  an  interactive  program  for  use  at  a  time-sharing  demand 
terminal.  It  provides  the  user  with  the  equivalent  of  a  calculator  designed 
to  perform  signal  processing  operations  and  provides  the  software  equivalent 
of  a  large  number  of  modules  for  breadboarding  a  complete  signal  processing 
system. 

This  interpreter  is  meant  to  be  used  by  engineers  and  scientists  who  are 
familiar  with  signal  processing,  but  who  may  have  no  knowledge  of  programming. 

The  signal  processing  interpreter  is  particularly  useful  for  the  rapid 
investigation  of  systems  whose  complexity  precludes  a  complete  analytic  study, 
and  whose  utilization  of  new  components  may  make  hardware  breadboarding  un¬ 
desirable  because  of  cost  and  procurement  time  limitations. 

The  present  note  describes  software  additions  and  modifications  generated 
during  the  period  July  1973  through  March  1974  and  is  meant  to  be  used  in 
conjunction  with  the  Signal  Processing  Interpreter  Preliminary  Description 
and  User's  Guide  [1],  It  replaces  NUC  TN  1213,  by  the  same  name  [3], 
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A.  New  Commands  and  Newly  Implemented  Commands 

As  the  Signal  Processing  Interpreter  has  been  used  and  tested  by  potential 
and  actual  users,  the  need  for  extensions  and  expansions  of  the  system  have 
become  evident.  These  needs  were  mostly  met  by  the  addition  of  a  variety 
of  new  commands  to  the  repertory  of  existing  commands. 

Below  is  an  alphabetic  listing  of  these  new  additions  giving  the  names 
only.  Their  functions  and  use  will  be  covered  in  detail  in  the  body  of  this 
report,  under  Section  B. 

'ADD' 

'AUTOSCALE' 

'BARKER  GENERATOR1 
'CIRCULANT  PROJECTION' 

'CLEAR' 

'CLEAR  X' 

'COMPLEMENTARY  SEQUENCE' 

'CORRECT' 

'CORRECT  MACRO' 

'DECIMATE' 

'DEFINE  DPCM  GAIN' 

'DEFINE  DPCM  QUANTIZER' 

'DEFINE  DPCM  RCVR' 

'DEFINE  MACRO' 

'  DFT' 

'DPCM' 

'DPCM  RCVR' 

' EDCT ' 

'EXECUTE  MACRO' 

'EXPAND' 

'EXPONENTIAL' 

' IMAG  PART' 

'LIST  DPCM' 

'LIST  MACRO' 

'MASS  RECALL' 

'MASS  STORE' 
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'MODIFIED  LOG' 

'NOISE  GENERATOR' 

' ODCT ' 

'PLOT' 

'POLAR  TO  RECTANGULAR' 
'PRINT  ALL  COMMANDS' 
'QUANTIZER' 

'REAL  PART' 

'RECIPROCAL' 
'RECTANGULAR  TO  POLAR' 
'RIGHT  CIRCULAR  SHIFT' 
'RIGHT  SHIFT' 

'ROLL  DOWN' 

'ROLL  UP' 

'SCALE' 

'SUM' 
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B.  Detailed  Description  of  the  New  Commands  and  Newly  Implemented  Commands* 

1 .  Exponentiation 

Exponentiation  of  the  elements  of  a  vector  was  made  available  to  users 
of  the  Signal  Processing  Interpreter  by  the  new  command 

'EXPONENTIAL' 

It  simply  replaces  each  component  Zi  of  the  vector  by  e^ ,  where  Zi 
is  complex. 

2.  Reciprocals 

The  comnand  'RECIPROCAL'  was  added.  It  checks  for  zero  data  points 
of  a  vector  to  avoid  illegal  division  formulations.  If  one  is  en¬ 
countered,  an  error  message  results  that  says: 

"VECTOR  CONTAINS  A  ZERO.  RECIPROCAL  NOT  PERFORMED." 

and,  in  fact,  no  reciprocal  is  then  performed.  If  all  data  points  are 
legitimately  non-zero,  then  processor  replaces  each  point  of  the  vector 
by  its  complex  reciprocal. 

3.  Addition  of  Vectors 

Addition  of  vectors  was  made  possible  by  including  the  command  'ADD'. 
This  will  cause  addition  of  the  lowest  two  vectors  in  the  stack,  after 
zero  filling  the  shorter  one  to  the  size  of  the  longer.  It  will  also 
lower  the  remaining  vectors  of  the  stack  one  level.  Resultant  is 
placed  in  the  X-register. 


*In  the  following  naragraphs  that  describe  the  commands,  use  is  made  of  a  set  of 
sub  paragraphs  under  the  headings  "Programming  Considerations".  These  were 
designed  for  the  person  who  must  know  in  detail  how  and  why  the  code  was  written 
the  way  it  was.  The  average  user  can  either  skim  them  or  skip  them  entirely. 
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4.  Shifts 


a.  Commands 

There  are  two  basic  kinds  of  shifts  available  to  the  user  of  the  Signal 
Processing  Interpreter.  There  is  a  straight  shift  of  the  elements  of 
a  vector,  with  the  end  values  being  lost  as  they  are  shifted  off  the  end 
and  zeros  filling  the  new  positions  opened  up.  This  is  called  by: 

'RIGHT  SHIFT' 

And  there  is  a  circular  shift  of  the  vector  components,  with  all  values 
being  preserved.  This  is  called  by: 

■RIGHT  CIRCULAR  SHIFT' 

There  are  no  'LEFT  SHIFT'  or  'LEFT  CIRCULAR  SHIFT'  commands.  These 
cases  are  handled  in  another  way.  (See  below). 

b.  Use 

(1)  To  perform  a  right  shift,  user  puts  his  vector  into  the  X-register 
and  calls  'RIGHT  SHIFT'.  Processor  then  interrogates  him  with  message 
"ENTER  NUMBER  OF  POSITIONS  TO  BE  SHIFTED",  to  wlich  he  responds  by 
entering  an  integer.  Vector  is  then  shifted  that  number  of  places 

to  the  right, while  that  number  of  vector  components  are  lost  off  the 
right  end.  Also  the  left  end  will  now  contain  that  number  of  zero 
components.  Shifted  vector  will  appear  in  the  X-register. 

(2)  To  perform  a  left  shift,  everything  is  done  the  same  way  except  that 
the  number  entered  by  the  user  for  positions  to  be  shifted,  is  entered 
as  a  negative  integer.  In  this  case  the  zero  components  appear  at  the 
right  end  and  shifted  components  are  lost  off  the  left  end. 

(3)  To  perform  a  right  circular  shift,  user  puts  his  vector  in  the  X-register 
and  calls  'RIGHT  CIRCULAR  SHIFT'.  He  receives  the  same  message 

"ENTER  NUMBER  OF  POSITIONS  TO  BE  SHIFTED" 
to  which  he  responds  by  entering  an  integer.  A  circular  shift  then 
takes  place  and  resultant  appears  in  the  X-register. 
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(4)  To  perform  a  left  circular  shift,  everything  is  done  the  same  way  except 
that  the  number  entered  by  the  user  for  positions  to  be  shifted  is 
entered  as  a  negative  integer. 

(5)  If  the  number  entered  for  either  shift  command  is  zero,  nothing  takes 
place  and  user  is  asked  for  next  command. 

(6)  If  the  number  of  places  shifted  is  greater  than  or  equal  to  the  vector 
size  for  a  straight  shift,  the  result  is  to  zero  out  your  vector, 

.  Programming  Considerations 

(1)  Circular  shifts  are  done  modulo  the  vector  size.  A  circular  shift 
by  the  vector  size  leaves  the  vector  unchanged. 

(2)  Circular  shifts  are  done  by  a  FORTRAN  subroutine  called  CSHIIT,  while 
straight  shifts  are  done  by  Code  embedded  in  ALGOL  routine  SPIN  3. 

(3)  In  order  not  to  cause  problems,  the  value  used  for  a  straight  shift 
is  made  equal  to  the  smaller  of  the  following  two  values:  vector  size 


and  number  entered. 

(4)  Since  a  circular  shift  by  a  multiple  of  NX  is  equivalent  to  a  shift  of 
zero,  a  circular  shift  of  -k  is  equivalent  to  a  circular  shift  of  NX-k, 
so  the  distinction  between  right  circular  shifts  and  left  circular  shifts 
is  purely  conceptual. 


5-  Stack  Operations 

a.  Commands 

Several  new  stack  operations  were  provided  arid  are  listed  below: 

'ROLL  UP* 

'ROLL  DOWN' 

'CLEAR' 

'CLEAR  X' 

b.  Us£ 

(1)  User  would  use  these  in  same  way  as  the  stack  operations  on  the  HP-35. 

Their  usage  is  straightforward.  See  reference  [7]  for  information  on  the  HP-35. 

(2)  'CLEAR'  initializes  with  zero  all  512  cells  of  all  registers^ regardless 
of  their  original  vector  dimensions.  Then  it  sets  all  register  dimen¬ 
sions  to  the  nominal  value  of  1.  'CLEAR  X'  does  so  only  for  the  X-register. 

(3)  'ROLL  UP'  and  'ROLL  DOWN'  also  work  on  512  cell  transfers  regardless 
of  vector  size. 

c.  Programming  Considerations 

(1)  The  vectors  were  cleared  as  512  and  rolled  up  or  down  as  512  for  two 
reasons.  One  was  efficiency  of  code  written  and  the  other  to  avoid 
leaving  garbage  in  cells  that  might  give  trouble  at  a  later  time. 

(2)  'ROLL  UP'  and  'ROLL  DOWN'  actually  perform  cyclic  permutations  of  the 
vectors  X,  Y,  Z,  T  and  also  of  their  dimensions  NX,  NY,  NZ  and  NT. 
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6.  Seal i nq 

a.  Commands 


Scaling  can  be  done  to  a  vector  by  using  one  of  the  following  commands: 
'SCALE' 

' AUTOSCALE' 

The  first  one  requires  a  scale  factor  to  be  supplied  by  the  user,  while 
the  second  one  is  an  automatic  scaling  process. 

b.  Use 

(1)  The  commard  'SCALE'  calls  upon  the  user  to 

'ENTER  COMPLEX  SCALE  FACTOR' 

When  he  does,  his  original  vector  is  multiplied  pointwise  by  this  scale 
factor,  and  the  result  is  in  the  X-register. 

(2)  When  'AUTOSCALE'  is  called,  processor  calculates  the  absolute  value 
of  the  largest  point  of  the  original  vector  and  scales  vector  by  this 
quantity.  It  also  prints  out  the  message 

"PREVIOUS  LARGEST  ABSOLUTE  VALUE  =  _ « 

The  variable  printed  is  this  quantity  used  for  scaling, 

c.  Programming  Considerations 

(1)  A  FORTRAN  subroutine  called  SCALE  was  used  even  though  the  code  was  so 
short  to  suppose  that  it  would  logically  be  embedded  in  the  ALGOL  SPIN  3 
procedure.  This  was  done  to  allow  uniform  entry  of  the  complex  scaling 
factor.  To  enter  complex  numbers  by  FORTRAN  is  different  from  ALGOL. 

An  entry  of  a  complex  number  for  FORTRAN  is  of  the  form 

4.,  -3.6 

while  for  ALGOL,  it  is  of  the  form 

4.  -3.6 

(2)  Autoscaling  makes  use  of  the  ALGOL  procedure  MAX. 


•  Correction  of  Vector  Points 

a.  Command 


The  user  may  correct  a  point  of  a  vector  in  the  X- 
' CORRECT' 


register  by  calling 


b.  Use 


0)  After  user  places  his  call  to  'CORRECT',  processor  comes  back  to  ask 
user  to 


"ENTER  INDEX  OF  POINT  TO  BE  CORRECTED" 

(2)  User  then  enters  his  integer  and  receives  back  the  message 

"ENTER  COMPLEX  DATA  POINT  NUMBER 

- -  > 

where  the  variable  is  the  index  he  just  entered. 

(3)  Now  he  enters  a  complex  value  in  FORTRAN  format  and  receives  a  verified 
tion  print  back  of  what  he  just  entered,  with  the  message 

"VALUE  STORED  AT  INDEX  IS 

■  -  -  _ _ • 

(4)  Processor  replaces  previous  value  at  that  point  by  the  new  one. 

c.  Programming  Considerations 


(1)  Correction  is  done  in  a  FORTRAN  subroutine  named  CORRV,  This  was  done 
to  allow  use  of  standard  FORTRAN  format  for  inputing  the  complex  value 
of  the  point  to  be  corrected.  ALGOL  uses  a  different  format,  and  it  was 
necessary  to  have  consistency  of  input  conventions. 

(2)  At  some  later  time.  It  may  be  desirable  to  make  sure  index  entered  is 

within  the  bounds  of  the  vector  being  corrected.  No  check  is  now  made 
and  problems  could  result  from  sloppy  use. 


(3)  complex  input  format  for  FORTRAN  has  the  real  part  first,  then  a  comma, 
then  the  imaginary  part.  The  real  part  and  the  imaginary  parts  are 
rea 1  numbers  with  decimal  points. 
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8.  Printing  Command  Names 

a.  Command 

An  operator  was  implemented  to  allow  the  printing  out  of  the  command 
name  as  soon  as  it  is  recognized.  It  is  called  by 

'PRINT  ALL  COMMANDS'  . 

b.  Use 

(1)  Whenever,  if  ever,  the  user  desires  to  see  the  command  names  printed 
out,  he  need  only  call  once,  'PRINT  ALL  COMMANDS',  and  it  will  be 
done  for  the  remainder  of  the  run. 

(2)  This  was  implemented  mainly  for  batch  jobs  where  a  sequential  listing 
of  the  commands  performed  is  otherwise  unavailable. 

(3)  This  printing  of  commands  will  also  occur  for  the  elements  that  make 
up  a  macro. 

c.  Programming  Considerations 

(1)  Processor  uses  a  flag  called  PCFLAG.  It  has  an  initial  value  of  0 
and  changes  to  1  at  the  first  time  'PRINT  ALL  COMMANDS'  is  called. 
Thereafter  it  remains  1  for  the  rest  of  the  run.  In  other  words, 
there  is  no  way  to  shut  it  off,  at  present.  An  additional  command 
could  easily  be  created  to  do  that.  It  need  only  reset  the  flag. 

(2)  This  command  provides  a  printout  for  batch  jobs  which  is  readable 
without  reference  to  the  card  deck  which  was  used  for  submission  of 
the  job. 


9.  Polar  and  Rectangular  Coordinate  Conversion 

a.  Commands 


User  can  perform  conversion  of  a  vector  from  polar  to  rectangular 
coordinates  by  calling: 

'POLAR  TO  RECTANGULAR' 

Or  he  can  convert  rectangular  coordinates  to  polar  by: 

'RECTANGULAR  TO  POLAR' 

b.  Use 

When  either  of  the  above  is  called,  the  user  receives  message: 

"ENTER  R  OR  D  IN  STRING  QUOTES  FOR  ANGLES  IN  DEG,  RAD" 

User  responds  with  ' D'  or  ' R '  . 

c.  Programming  Considerations 

(1)  Vector  data  is  complex  type,  so  when  rectangular,  the  real  part  is 
X  coordinate  and  imaginary  part  is  Y  coordinate.  For  polar,  real 
part  is  magnitude  and  the  imaginary  part  is  angle. 

(2)  Program  actually  only  checks  for  the  ' D '  option.  If  not  ’ D ' ,  then 
' R '  is  assumed. 

(3)  It  was  necessary  to  include  a  special  routine  to  compute  the  angles. 

The  ALGOL  arctangent  routine  does  not  preserve  uniqueness  through  the 
transform  since  it  uses  only  a  single  input  argument.  The  FORTRAN 
routine  FUNCTION  CANG  was  added.  It  also  handles  the  case  of  zero 

for  both  coordinates.  To  use  it  however,  it  was  necessary  to  type  it  as 

EXTERNAL  FORTRAN  REAL  PROCEDURE  . 
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a  .  Command 


User  may  call  'OFT'  when  a  discrete  Fourier  transform  is  desired  and  not 
have  to  make  the  decision  himself  about  whether  it  is  an  * FFT 1  or  a  ' CZT ' . 

b.  Use 

User  puts  vector  to  be  transformed  into  the  X-register,  calls  'DFT',  and 
receives  his  transformed  vector  back  in  the  X-register. 

c.  Programming  Considerations 

(1)  Program  will  check  the  vector  dimension  for  powers  of  2  up  to  an  exponent 
of  9,  that  is,  up  to  a  value  of  512.  If  a  power  of  2,  then  the  FFTSP 
subroutine  will  be  called.  Otherwise  it  will  be  CZT. 

(2)  Code  is  redundant  with  some  in  ' EDCT ' .  (See  writeup  of  Cosine  Trans¬ 
forms).  It  may  be  desirable  to  convert  this  to  a  subroutine  and  replace 
the  redundant  code  with  subroutine  calls. 


1 1. Noise  Generation 

a.  Command 

A  noise  generation  capability  is  available  to  the  Signal  Processing 
Interpreter  user  by  calling 

'NOISE  GENERATOR' 

b.  Use 

(1)  After  calling  'NOISE  GENERATOR',  user  receives  back  the  message: 

"ENTER  INITIALIZATION  POINT  FOR  NOISE  GENERATOR"  . 

(2)  He  responds  by  entering  a  positive  integer  that  controls  the  starting 
point  of  the  random  number  generator  that  creates  the  noise  sequence. 

(3)  User  then  is  asked  to:  "ENTER  VARIANCE'’.  This  is  a  real  number  and 
represents  the  variance  that  the  final  noise  sequence  will  have. 

(4)  After  he  enters  this  number,  the  noise  sequence  is  generated  and 
appears  in  the  X-register. 

(5)  The  length  of  the  noise  sequence  generated  is  the  length  of  the  X-register 
vector  at  the  time  it  is  called. 

c.  Programming  Considerations 

(1)  Processor  calls  upon  a  FORTRAN  subroutine  called  NOISE  to  perform  the 
task. 

(2)  The  method  used  for  the  random  number  generator  came  from  "Communications 
of  the  ACM",  Volume  12,  Number  2,  of  February  1969,  p.  93-94,  by 

J.  B.  Kruskal.  [2] 

(3)  The  random  number  generator  used  was  set  up  to  allow  sequences  of  up 
to  8192  points.  Since  we  are  limited  to  512,  in  general,  we  can  start 
the  sequences  at  a  variety  of  points  to  get  independent  sequences. 

Method  generates  a  zero-mean  random  vector  using  a  multiplicative  con- 
gruential  method. 
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a.  Commands 

Expansion  and  decimation  of  a  vector  are  available  to  the  user  through 
the  commands: 


'EXPAND' 

'DECIMATE' 


b.  Use 

(1)  To  expand  a  vector  in  the  X-register,  user  calls  'EXPAND'.  Processor 
then  asks  user  to 


"ENTER  INTEGER  EXPANSION  FACTOR"  , 


to  which  the  user  responds  with  a  positive  integer.  (Call  it  INT) 

(2)  The  new  vector  will  be  of  length 

1  +  ((NX-1)  *  INT),  where  NX  was  the  original  length.  The  first 
point  remains  fixed,  while  between  each  two  points  of  the  original 
vector  there  will  be  inserted  (INT-1)  zero  points.  Resultant  appears 
in  the  X-register  and  destroys  the  original  vector. 

(3)  Decimation  works  in  a  corresponding  way.  With  the  vector  in  the  X-regijter, 
the  user  calls  'DECIMATE'.  The  processor  then  asks  user  to 

"ENTER  INTEGER  DECIMATION  FACTOR"  , 
to  which  the  user  responds  with  a  positive  integer.  (Call  it  NNT). 

(4)  The  new  vector  will  be  of  length  1  +  ( ( NX-1 )/NNT )  where  MX  was  the 
original  vector  size.  The  new  vector  will  be  the  first  point  and  every 
NNT-th  point  thereafter. 

(5)  Note  that  these  operators  are  not  commutative,  either  in  value  or 
in  length. 


Programming  Considerations 


(1)  The  code  for  decimation  uses  integer  division  as  noted  by  the  //  sign. 
This  is  necessary  to  truncate  off  fractional  parts. 
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13.  Cosine  Transforms 

a.  Commands 

Even  and  odd  discrete  cosine  transforms  are  available  through  the  commands 

' EDCT ' 

'  ODCT ' 

b.  Us£ 

(1)  User  puts  vector  to  be  transformed  in  the  X-register  and  then  gives 
appropriate  call:  'EDCT'  for  even  discrete  cosine  transform  or  'ODCT' 
for  odd  discrete  cosine  transform. 

(2)  Processor  generates  transform  and  stores  the  results  back  into  the 
X-register. 

(3)  Size  of  the  original  vector  is  limited  to  256,  with  the  one  exception 
of  the  number  512  being  allowed  for'the  'EDCT'. 

c.  Programming  Considerations 

(1)  Actual  formulas  implemented  are  not  included  here  nor  is  an  explanation 
of  the  method.  User  is  referred  to  NUC  TN  1265,  entitled  "High  Speed 
Serial  Access  Implementation  for  Discrete  Cosine  Transforms"  by  J.  M. 
Speiser. 

(2)  The  reason  for  limitation  of  size  of  original  vector  is  due  to  method 
used.  If  not  a  power  of  2,  method  makes  use  of  CZT  routine  which 
quadruples  size.  But  the  dimensioning  is  limited  to  1024,  thus  256 
for  original  vector. 

(3)  Routines  make  use  of  a  flag  called  NODD, which  takes  on  the  value  1  for 
odd  or  0  for  even.  The  odd  transform  automatically  uses  the  CZT  routine, 
while  the  even  one  must  test  for  power  of  2.  It  it  is  a  power  of  2, 
then  the  FFT  routine  is  used;  otherwise  the  CZT  is  again  used. 

(4)  We  cannot  use  the  ' DFT '  formal  call  to  do  this  as  the  code  is  embedded 
in  the  SPIN  3  routine,  but  the  code  to  do  this  checking  is  identical. 
One,  at  some  later  time,  could  make  the  'DFT'  a  subroutine  and  then 
use  it  for  both  places^as  the  code  is  redundant. 
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14.  Modified  Log 

a.  Command 

A  modified  logarithm  routine  is  available  to  user  with  call 
’MODIFIED  LOG' 

b.  Use 

(1)  Vector  to  be  processed  is  placed  in  the  X-register  and  call  is  made 
to  'MODIFIED  LOG' 

(2)  At  this  point  program  interrogates  user  with  message 

"ENTER  MINIMUM  RECOGNIZED  NUMBER' 

(3)  User  then  enters  a  positive,  real  number.  If  it  is  not  positive, 
a  diagnostic  results 

"ERROR.  NUMBER  MUST  BE  POSITIVE  REAL" 
and  he  is  again  asked  to  enter  his  number  until  he  gets  it  right. 

(4)  The  real  part  of  the  vector  is  then  checked  against  this  minimum 
recognized  number.  If  it  is  larger,  then  the  complex  number  is  replaced 
by  the  natural  logarithm  of  its  real  part.  If  not,  it  is  replaced  by 
the  natural  logarithm  of  the  minimum  recognized  number. 

c.  Programming  Considerations 

(1)  Code  will  not  at  present  catch  errors  in  type  for  entering  the  minimum 
recognized  number.  But  this  safeguard  could  easily  be  built  in. 
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15.  MACROS 


a.  Commands 

User  may  define  and  use  combinations  of  Signal  Processing  Interpreter 
statements  as  blocks  or  macros.  To  do  so  he  will  use  some  or  all  of  the 
following  four  commands: 

'DEFINE  MACRO' 

'LIST  MACRO' 

'CORRECT  MACRO' 

'EXECUTE  MACRO' 

b.  Use 

(1)  To  begin  the  process,  user  enters  'DEFINE  MACRO'  .  Only  one  macro  may  be  used 
at  one  time  by  the  user.  The  processor  then  prints  out: 

"ENTER  NUMBER  OF  COMMANDS  TO  BE  USED  FOR  MACRO." 

The  user  then  enters  an  integer  between  2  and  20.  (Actually  1  would 
work,  but  would  be  meaningless).  If  he  enters  a  number  greater  than 
20,  a  diagnostic  results. 

"ERROR.  MORE  THAN  20  COMMANDS" 

He  is  again  asked  for  number  of  commands. 

(2)  Once  he  has  properly  entered  the  number,  he  receives  the  message: 

"ENTER  COMMAND"  until  he  has  entered  enough  statements  to  make 
up  the  desired  size  for  the  macro.  At  this  point  he  receives  message: 

"END  OF  MACRO  DEFINITION" 

(3)  Provisions  are  made  to  keep  the  run  from  terminating  if  user  makes  a 
bad  entry.  He  then  receives  a  message: 

"ERROR  IN  ENTERING  LAST  MACRO  COMMAND" 
and  is  given  another  chance  to  enter  it  correctly. 

Note:  Processor  does  not  do  a  command  table  search  at  this  time  so  command 
still  may  not  be  legitimate.  It  just  means  that  entry  was  of  correct 
format. 
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(4)  User  may  check  to  see  what  Is  in  the  macro  by  calling: 

•LIST  MACRO' 

It  gives  a  printed  listing  of  the  command  names  that  make  up  the  macro, 
and  tnis  process  makes  a  good  check. 

(5)  The  user  may  desire  to  correct  the  generated  macro  or  change  it.  He 
may  do  so  by  calling 

'CORRECT  MACRO' 

This  is  similar  in  concept  to  the  command  'CORRECT'  associated  with  changes 

in  vector  values. 

Program  returns  message: 

"ENTER  INDEX  OF  COMMAND  TO  BE  CORRECTED" 

User  enters  an  integer  between  1  and  20  and  receives  back  message: 

“ENTER  COMMAND"  ,  which  he  does.  Process  repeats  if  an  error  occurs 
in  entering.  Upon  successful  completion,  control  returns  to  normal 

sequence. 

(6)  User  is  now  ready  to  execute  his  macro,  done  by  calling: 

•EXECUTE  MACRO' 

This  will  cause  the  first  command  to  be  executed  and  all  subsequent 
others  possible  up  to  the  first  required  parameter  input  or  option 
selection.  Process  continues  until  all  commands  in  macro  have  been 
executed,  at  which  time  user  receives  message: 

"END  OF  MACRO  EXECUTION" 

(7)  Control  returns  to  normal  sequence. 

(8)  User  note  that  the  normal  sequence  message  of  the  form 

"NX=8  PLEASE  ENTER  NEXT  COMMAND" 

will  be  deleted  between  executed  commands  of  a  macro.  It  will  aopear 
after  the  completion  of  any  one  of  the  four  macro  control  commands. 

(9)  Macro,  once  generated,  is  available  to  the  user,  at  his  discretion, 
for  the  rest  of  the  run. 
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Programming  Considerations 

(1)  Macros  are  limited  to  20  commands  at  present.  This  could  easily  be 
changed  by  minimal  reprogramming.  Also  more  than  one  macro  could  be 
implemented  by  a  somewhat  larger  job  of  programming. 

(2)  Read  error  checks  are  built  in  to  catch  errors  in  entering  commands 

to  macro.  Check  is  made  for  form  only.  Errors  in  spelling  or  erroneous 
names  won't  be  caught  until  user  tries  to  execute  his  macro. 

(3)  A  flag  named  MAFLAG  is  used  to  control  suppression  of  unwanted  print 
statements  between  macro  commands.  It  is  set  to  1  during  execution 
of  a  macro  and  to  0  at  all  other  times. 


a.  Commands 


Mass  memory  for  vector  storage  and  retrieval  is  available  to  the  signal 
processing  interpreter  user  through  the  following  commands: 

'MASS  STORE' 

'MASS  RECALL' 

b.  Use 

(1)  In  order  to  use  this  storage  capability,  the  user  must  include  in 
his  preliminary  setup  commands  (before  execution)  the  following  two 
statements  (or  variations  of  them): 

@ASG,  AZ  NAME. 

OUSE  15,  NAME. 

where  NAME  is  the  users  own  catalogued  file,  previously  assigned  by: 

OASG,  UP  NAME.,  F/120/TRK 

The  120  tracks  of  mass  storage  file  are  necessary  to  allow  for  a  maximum 
of  201  vectors  of  512  complex  words.  User  should  supply  his  own  file 
name  for  NAME. 

(2)  Now  the  Signal  Processing  Interpreter  can  be  executed. 

(3)  The  two  storage  commands  are  very  similar  in  use.  To  store  a  vector, 
the  user  puts  his  vector  into  the  X-register  and  types  'MASS  STORE', 

The  Signal  Processing  Interpreter  then  interrogates  him  with  the 
message: 

"ENTER  INDEX  OF  MASS-STORED  VECTOR"  . 

User  then  enters  an  integer  between  0  and  200.  The  program  will  store 
the  vector  in  the  appropriate  memory  area  along  with  its  size.  Con¬ 
versely,  to  retrieve  a  vector,  user  types 

'MASS  RECALL' 

at  which  time  he  receives  the  same  interrogation  as  above.  When  he  enters 

the  integer  (0  <  N  <  200),  he  receives  the  vector  back  in  the  X-register 

along  with  its  size  in  NX.  User  can  store  up  to  201  vectors  at  a  time. 

'MASS  RECALL'  treats  the  stack  in  the  same  way  as  'RECALL',  that  is,  it 

raises  the  stack  one  level  with  the  topmost  vector  being  lost. 
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c.  Programming  Considerations 


(1)  Program  uses  same  code  to  accomplish  both  commands, with  a  read-write 
flag  as  the  only  difference. 

RWFLAG  =  1  for  write  and  RWFLAG=2  for  read. 

(2)  Program  makes  access  to  mass  storage  by  sector  addressing.  In  order 
to  do  so,  program  computes  and  locates  proper  spot  by  using  internal 
calculations  in  SUBROUTINE  MASTOR  and  by  calling  system  SUBROUTINE 
SETADR  with  correct  index. 

(3)  To  make  the  above  calculation  one  needs  to  know  that  there  are  28  words 
per  sector.  In  order  to  be  able  to  store  a  512  complex  word  vector, 

we  would  need  at  least  1024  words.  This  will  be  handled  by  37  sectors 
(1036  words).  In  addition,  the  program  writes  and  stores  a  one-word 
block,  before  storing  the  vector,  that  gives  vector  length.  This  reauires 
an  additional  sector  for  a  total  of  38.* 

(4)  If  an  index  of  200  (actually  the  201st  position)  is  exceeded,  user  is 
given  a  diagnostic  message: 

"EXCEEDED  MAX  INDEX.  OPERATION  NOT  DONE." 

Then  the  mass  storage  operation  is  not  performed,  with  control  returning 
to  ask  the  user  for  the  next  command. 

(5)  The  NTRAN  no-op  statement  before  CALL  SETADR  of  SUBROUTINE  MASTOR  is 
required  in  order  to  assign  unit  15  as  an  NTRAN  unit.  If  it  is  not 
present,  the  call  to  SETADR  assumes  15  is  a  FORTRAN  file,  and  then  the 
first  NTRAN  call  to  unit  15  causes  the  run  to  bomb.  (A  file  cannot  be 
both  FORTRAN  and  NTRAN). 


*There  are  64  sectors  to  1  track.  Since  we  need  201  vectors  each  of  38  sectors 
length,  we  need  a  total  of  7638  sectors,  or  about  120  tracks. 
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17.  Circulant  Projection 

a.  Command 

User  may  generate  a  circulant  projection  of  a  given  vector  by  calling: 
'CIRCULANT  PROJECTION' 

b.  Use 

(1)  This  replaces  Xp  by  (1  -  £)Xp  +  J  X*_p 

for  P  =  0,  1,  2,  N-l  and  where  *  means  complex  conjugate. 

(2)  User  puts  Hs  vector  to  be  processed  in  the  X-register,  calls  'CIRCULANT 
PROJECTION',  and  then  receives  this  projection  back  in  the  X-register. 

c.  Programming  Considerations 

(1)  Code  makes  use  of  several  procedures  already  defined,  as  well  as  several 
temporary  arrays  for  intermediate  results. 

(2)  Formula  above  supplied  by  Jeff  Speiser. 


18.  Barker  Sequence  Generation 

a.  Command 

User  may  generate  a  Barker  sequence  by  a  call 
'BARKER  SEQUENCE' 

b.  Use 

(1)  By  such  a  call,  a  Barker  sequence  will  be  generated  in  vector  form 
in  the  X-register. 

(2)  Processor  asks  user  to 

"SELECT  SIZE  FROM  2,  3,  4,  5,  7 ,  1 1 ,  1 3"  . 

(3)  User  then  enters  size  as  an  integer,  one  of  the  above.  If  he  makes  a 
mistake  and  enters  some  other  integer,  he  will  get  either  the  7  sequence 
or  the  13  sequence,  depending  upon  the  size  of  the  number  he  entered. 

But  he  is  not  kicked  off. 

c.  Programming  Considerations 

(1)  Barker  sequences  are  from  a  book  by  C.E.  Cook  and  M.  Bernfeld.  See 
reference  [5].  The  above  numbers  represent  all  known  sequences. 

(2)  Where  more  than  one  sequence  is  known  (i.e.  for  2  and  4),  only  one  of 
them  is  used. 

(3)  Code  for  generating  sequences  is  contained  in  a  FORTRAN  subroutine 
called  BARKER. 

(4)  In  order  to  conserve  core  storage,  sequences  are  overlapped  where 
possible.  The  variable  J  gives  the  starting  point  index  of  the  sequence 
in  the  input  data  array. 

(5)  Sequences  are  given  using  +1  and  -1  rather  than  just  +  and  -  as  in  paper. 
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19.  Complementary  Sequence  Generation 

a.  Command 

A  complementary  sequence  generation  capability  is  available  to  the  user 

by  a  call  to 

'COMPLEMENTARY  SEQUENCE' 

b.  Use 

(1)  User  can  generate  a  complementary  sequence  of  length  10,  26,  or  a 
power  of  2  up  to  9,  that  is,  the  number  up  512. 

(2)  After  calling  'COMPLEMENTARY  SEQUENCE',  the  user  receives  back  the 
message: 

"SELECT  SIZE  10,  26,  or  0-9  (FOR  POWER  OF  2)" 

(3)  If  he  selects  10  or  26,  a  pre-stored  sequence  will  be  fed  back  into 
the  X-register  with  the  appropriate  dimension.  If  it  is  a  power  of  2, 
then  the  sequence  is  generated  by  block  recursion:  A'  =  AR,  B  =  AB. 
See  reference  [6]. 

(4)  Sequences  form  vectors  of  complex  components, 

c.  Programming  Considerations 

(1)  Actual  sequence  generation  takes  place  in  a  FORTRAN  subroutine  called 
COMPLM. 

(2)  If  user  makes  a  mistake  and  chooses  some  other  integer  than  those 
allowed,  he  will  get  either  the  26  sequence  for  numbers  greater  than 
10  or  he  will  get  the  sequence  for  2  to  the  power  1  for  neoative 
values.  This  is  due  to  the  way  the  tests  are  made.  These  will  not 
be  what  he  asked  for.  but  he  will  not  have  his  run  terminated  either. 

(3)  Real  and  imaginary  parts  of  the  10  and  26  sequences  are  stored  in  real 
arrays  and  assembled  by  using  the  FORTRAN  function  CMPLX. 

(4)  Sequence  elements  are  made  up  of  +1  and  -1  rather  than  1  and  0  as 
given  in  reference  [6]. 
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20.  Differential  Pulse  Code  Modulation 


a.  Commands 

The  Signal  Processing  Interpreter  has  four  commands  which  make  possible 
the  definition  and  use  of  a  differential  pulse  code  modulation  and  two 
more  commands  for  an  associated  receive".  These  are: 

'DEFINE  DPCM  GAIN' 

'DEFINE  DPCM  QUANTIZER' 

'LIST  DPCM' 

'DPCM' 

'DEFINE  DPCM  RCVR' 

'DPCM  RCVR' 

b.  Usj; 

(1)  In  order  to  use  the  'DPCM'  command  or  the  'DPCM  RCVR'  command,  the 
user  must  first  define  then.  He  does  so  by  using  the  other  commands. 

If  no  definition  takes  place,  then  the  gain  is  assumed  to  be  of  magnitude  1 
and  the  quantizer,  if  present,  is  assumed  to  be  of  an  infinite  no.  of 
steps.  (That  is,  output  =  input.) 

(2)  User  calls  ‘DEFINE  DPCM  GAIN'  in  order  to  define  or  redefine  gain. 

Program  then  interrogates  user  with  message  "ENTER  LOOP  GAIN.  REAL 
NUMBER  BETWEEN  0.  AND  1."  ,  to  which  the  user  responds.  Note  that 
the  default  value  is  1. 

(3)  The  DPCM  is  still  considered  undefined  until  user  makes  first  call  to 
'DEFINE  DPCM  QUANTIZER'  .  Then  the  following  sequence  takes  place: 

(a)  Program  interrogates  user  with  message 

"DEFINE  QUANTIZER  1.  ENTER  NUMBER  OF  JUMP  POINTS" 

(b)  User  responds  with  an  integer  between  1  and  20. 

(c)  Program  then  asks  for  jump  points  to  be  entered  with  message: 

"ENTER  REAL  NUMBER  FOR  JUMP  _ "  . 

(d)  User  enters  jump  points  in  response. 

(e)  Same  sequence  as  (a),  (b),  (c),  and  (d)  for  Quantizer  2. 

(f)  DPCM  is  now  defined. 
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(4)  If  user  wants  to  see  how  DPCM  is  defined,  he  can  use 

'LIST  DPCM'  . 

This  will  print  out  gains  and  jump  points  for  the  two  quantizers. 

(5)  To  use  'DPCM  RCVR' ,  user  must  first  call  'DEFINE  DPCM  RCVR'  or  he  will 
get  the  default  case  of  gain  =  1. 

When  called,  it  interrogates  user  as  in  (2)  above.  When  answered, 

DPCM  RCVR  is  considered  defined. 

(6)  Once  defined,  user  need  only  call  'DPCM'  or  'DPCM  RCVR'  to  perform 
desired  operation. 

(7)  Redefinition  may  take  place  at  any  time  by  the  user. 

.  Programming  Considerations 

(1)  Internal  flags  are  set  the  first  time  either  DPCM  or  DPCM  RCVR  are 
defined.  After  that,  they  are  assumed  to  be  defined  for  the  rest  of 
the  run,  though  redefinition  may  take  place. 

(2)  Quantizers  are  limited  at  present  to  at  most  20  points.  But  more  could 
be  easily  added  with  simple  changes  to  the  program. 

(3)  When  user  tries  to  use  DPCM  or  DPCM  RCVR  but  has  not  defined  them,  he 
gets  a  diagnostic  message  along  with  having  program  use  default  conditions 

(4)  If  a  run  bombs  off  before  completion  and  execution  must  be  restarted, 
then  redefinition  must  take  place  again.  No  provision  is  made  for 
storing  DPCM  values  or  those  for  DPCM  RCVR  from  run  to  run.  Each  run 
reinitializes  all  quantities  involved. 


21 .  Quantization 


a.  Command 

User  has  available  a  quantizer  by  simply  calling 
'QUANTIZER'  . 

b.  Use 

(1)  User  puts  his  original  vector  in  the  X-register  and  calls  upon 
'QUANTIZER'.  Processor  immediately  asks  him  to: 

"ENTER  INTEGER  FOR  QUANTIZER" 

(2)  As  soon  as  user  does  so,  his  quantized  vector  is  generated  and  entered 
back  into  the  X-register. 

(3)  This  command  provides  a  uniform  quantization  into  the  specified  number 
of  magnitude  levels  with  preservation  of  sign.  Quantization  step  size 
is  equal  to  the  maximum  of  the  real  and  the  imaginary  parts  in  the 
array  divided  by  the  number  of  levels. 

c.  Programming  Considerations 

(1)  Processor  makes  use  of  a  FORTRAN  subroutine  called  OUANT  to  do 
the  quantization. 

(2)  The  scaling  factor  used  is  calculated  by  checking  all  the  components 
of  the  vector  for  the  largest  real  or  imaginary  component  magnitude. 
Then  the  integer  entered  is  divided  by  this  maximum  magnitude  to  make 
up  the  scaling  factor. 

(3)  Much  of  the  code  in  QUANT  involves  getting  the  rounding  process 
correct.  It  must  take  place  separately  for  real  and  imaginary  parts 
and  also  handle  negative  values.  The  real  and  imaginary  parts  are 
rounded  to  the  nearest  integer. 

22.  Real  and  Imaginary  Parts  of  Vectors 

Commands  to  extract  the  real  part  or  the  imaginary  part  of  a  complex 
vector  have  been  introduced  with  the  calls: 

'REAL  PART' 

' IMAG  PART' 

Note,  however,  that  'IMAG  PART'  replaces  the  X  vector  by  its  imaginary 
part  but  is  stored  then  as  a  real  vector.  That  is,  the  original  real  com¬ 
ponents  are  destroyed  and  the  original  imaginary  components  take  their 
place  by  moving. 
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a.  Command 


User  may  set  vector  dimension  by  using: 

'SET  X  DIMENSION' 

b.  Use 

(1)  After  giving  command  call,  user  receives  back  the  message. 

"ENTER  DESIRED  DIMENSION" 

To  answer  this,  he  enters  an  integer  between  1  and  512. 

(2)  If  a  nistake  is  made  entering  the  integer,  such  as  an  alphabetic 
character,  then  user  is  given  diagnostic  message. 

"ERROR  MADE  IN  ENTERING  DIMENSION.  TRY  AGAIN." 

And  sequence  starts  over  again  with  message  to: 

"ENTER  DESIRED  DIMENSION" 

(3)  When  dimension  is  properly  entered,  control  returns  to  ask  for  next 
command. 

c.  Programming  Considerations 

(1)  Code  makes  use  of  error  conditions  for  a  READ  statement.  These  will 
handle  almost  all  errors  except  an  entry  with  an  P  siqn,  signifying 

an  1108  control  statement.  This  causes  a  drop  in  level  and  termination 

of  execution. 

(2)  In  order  to  use  the  error  condition,  it  is  necessary  to  define  an  error 
exit  label.  One  cannot  just  define  it  with  a  colon  but  must  also  in¬ 
clude  the  label  name  among  the  non-executable  statements,  under  LOCAL 

LABEL. 


Inputting  Vectors  By  Points 

a.  Several  changes  were  made  to  facilitate  the  command  ’INPUT  VECTCR' . 

The  interrogation  of  the  user  was  transferred  entirely  to  the  FORTRAN 

subroutine  called  READV. 

b.  Next  READ  error  checks  were  built  into  the  reading-in  of  the  complex 
data  points.  Errors  caused  by  bad  inputting  of  data  will  no  longer 
terminate  a  computer  run  but  will  simply  cause  a  diagnostic  message 
to  be  printed  out. 

"ERROR  IN  ENTERING  DATA  POINT"  . 

Then  the  user  will  be  given  another  chance  to  enter  his  point  correctly. 
Most  errors  will  be  caught  but  not  ones  in  which  the  user  beqins  with 
a  master  space,  i.e.  0.  These  will  cause  a  drop  in  level,  and  user  must 
begin  his  run  all  over  again. 

c.  The  effect  of  this  is  that  noise  on  the  telephone  link  or  user  format 
errors  will  now  only  result  in  an  error  message  and  a  request  to  re¬ 
enter  the  command  or  data  point. 

d.  Programming  note:  it  became  necessary  to  read  the  subscript  of  the  input 
variable  with  the  variable  value  II+l  rather  than  the  normal  loop  variable 
value  I  which  you  would  expect.  This  was  due  to  a  FORTRAN  compiler  error 
discovered  while  trying  to  use  I.  Details  as  to  what  the  error  did  will 
not  be  covered  here.  It  suffices  to  say  that  the  error  was  verified 

by  the  System  Programming  Group  and  forwarded  on  to  UN  I VAC  for  correc¬ 
tion  via  an  SSFR  by  Ben  Ermance. 


3.  Plotting  Vectors 

There  are  now  two  recognized  calls  to  obtain  plots  of  the  X  vector.  Roth 
call  upon  the  same  code,  and  the  second  call  was  created  only  for  the 
convenience  of  the  user.  They  are: 

'PLOT  COMPLEX  VECTOR' 

'PLOT' 

See  reference  [1]  for  more  information. 

4.  Constant  Functions 

The  routine  'CONSTANT  FUNCTION'  was  rewritten  to  allow  a  consistent  con¬ 
vention  for  input  of  complex  values.  Delimiter  is  now  a  coma,  not  a 
blank  space  as  before.  If  more  details  are  needed,  see  paragraph  B.6.c.(l). 
The  constant  function  merely  fills  the  entire  vector  with  a  user-supplied 
complex  constant. 

5.  Printing  Vectors 

In  the  past  a  number  of  routines  automatically  produced  printed  listinas 
of  the  present  vector  contents.  This  has  been  changed  to  a  philosophy 
of  only  doing  what  the  user  asks.  Now  only  two  routines  will  result  in 
a  printed  listing: 

'PRINT  VECTOR' 

'PLOT'  or  'PLOT  COMPLEX  VECTOR' 

In  addition,  when  using  'PRINT  VECTOR',  you  no  longer  get  back  the 
message  saying  "FINISHED  PRINTING  VECTOR."  It  was  superfluous. 

6.  Other  Minor  Changes 

a.  Normalization  of  the  FFT  was  by  the  square  root  of  N.  The  CZT  uses  the 
FFTSP  subroutine,  and  thus  its  normalization  is  automatic.  Normalization 
is  not  now  done  within  the  CZT  routine  itself.  Due  to  possible  scaling 
problems,  this  procedure  will  probably  be  changed  in  the  near  future. 

b.  The  following  routines  were  rewritten  for  more  efficient  coding  but 
have  no  substantial  changes  to  their  basic  logic: 
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SUBROUTINE  CHIRP 
SUBROUTINE  CPLOT 
PROCEDURE  CZT 
SUBROUTINE  REV 
SUBROUTINE  SEE 
SUBROUTINE  TONE 

One  change  of  note  was  the  addition  of  EQUIVALENCE  statements  to 
SUBROUTINE  CPLOT  to  eliminate  the  need  for  the  FORTRAN  functions  REAL 
and  AIMAG.  This  is  far  more  efficient. 
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D.  Changes  to  SPIN3  Processor  System 


A  number  of  changes  were  made  to  the  SPIN3  procedure,  and  the  major  ones 
will  be  mentioned  below. 

a.  Previous  input/output  inconsistences  due  to  the  differences  between 
FORTRAN  and  ALGOL  conventions  have  been  eliminated.  This  mainly  affected 
the  inputting  of  complex  data  points. 

b.  All  vector  indices,  as  far  as  the  user  is  concerned,  are  row  consistent. 
They  run  from  zero  to  N-l  and  appear  in  inputting,  listing  and  plotting 
of  vectors. 

c.  Errors  that  previously  occurred  when  entering  commands  will  no  longer 
bomb  the  user  off  the  machine.  READ  error  checks  were  built  into  the 
system  to  handle  these  cases.  This  kind  of  error  now  causes  the  printing 
out  of  an  error  diagnostic  message 

'COMMAND  NOT  RECOGNIZED" 

and  then  a  return  to  ask  the  user  to  enter  his  next  command.  One  error 
that  still  is  not  caught  ar.d  corrected  is  when  user  includes  the  master 
space,  i.e.  @,  as  part  of  the  command.  This  will  cause  a  drop  in  level 
and  termination  of  the  run  being  executed. 

d.  The  output  printing  has  been  cleaned  up  to  minimize  unnecessary  line  feeds, 
inefficient  print  statements  and  the  printing  of  redundant  information. 
Controls  have  been  put  into  the  processor,  via  the  MARGIN  procedure  of 

the  ALGOL  processor,  to  delete  all  blank  lines  at  the  top  and  bottom  of 
the  print  pages.  This  was  done  to  allow  plots  to  be  conti nuous^ as  well 
as  listings.  This  change  affects  batch  jobs  only.  Clearly  cathode  ray 
scope  terminals  are  not  affected,  but  neither  are  the  DCT  500  terminals, 
unfortunately.  This  is  due  to  line-feeds  hard-wired  into  the  DCT  500 
hardware.  The  MARGIN  call  is 

MARGIN  ('M, 66, 0,0  .  ’) 

The  M  option  is  to  readjust  page  length  and  top  and  bottom  margins.  The 

0 

66  is  for  a  standard  page,  and  the  two  zeros  represents  lines  to  be  skipped 
at  top  and  bottom  of  page.  The  string  quotes  are  necessary  since,  in 
reality,  the  operator  MARGIN  has  only  one  argument,  namely  a  control 
string.  And  string  must  be  terminated  by  a  period.  This  takes  place 
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at  the  time  all  initialization  for  the  processor  is  done.  See  reference  [4] 
for  more  details  on  use  of  MARGIN. 

e.  Much  of  the  internal  code  has  been  rewritten  for  efficiency  of  execution 
and  compactness  of  format.  Thus  even  though  the  capabilities  of  the 
System  have  more  than  tripled,  the  number  of  punch  cards  has  remained 
roughly  the  same.  All  the  DO  loops  of  SPIN3  were  redone  in  compact  form 
as  were  many  IF  statements,  formats,  etc.  Procedures  internal  to  SPIN3, 
called  DOW,  PAD,  UP,  and  COPY  were  greatly  compacted. 

f.  An  initialization  process  was  introduced  into  the  processor  to  allow  a 
variety  of  operations  to  be  performed  only  once  during  the  run  to  eet  up 
the  conditions  and  prepare  for  execution  of  the  remaining  parts.  Prior 
to  this,  every  time  a  new  command  was  input,  the  processor  went  back  and 
reloaded  all  its  command  tables,  in  an  ALGOL  procedure  called  ATAB.  It 
then  went  through  an  elaborate  process  involving  Boolean  logical  operations 
to  perform  a  test  for  matching.  ATAB  has  been  eliminated  now,  the  command 
tables  are  loaded  only  once,  the  Boolean  operations  have  been  replaced  by 

a  simple  IF  statement,  and  all  initialization  is  now  done  only  once.  This 
was  done  by  putting  much  of  the  ATAB  code  into  the  SPIN3  processor  itself, 
by  making  all  initialization  the  first  executable  code  in  the  routine  and 
by  always  returning  to  points  below  this  code.  These  changes  improved 
running  times  for  the  executable  portion  of  the  program  by  a  factor  of 
between  4  and  5  times.  Another  change  contributing  to  this  great  time 
saving  was  terminating  the  table  search  as  soon  as  a  match  was  made.  Pre¬ 
viously  the  entire  table  was  searched,  even  if  recognition  had  already  taken 
place. 

g.  The  basic  flow  of  the  processing  of  a  command  is  now  slightly  more  complex 
since  the  case  of  macro  commands  must  be  considered.  Commands  that  make 
up  a  macro  take  a  different  path  upon  being  recognized. 

h.  Extensive  use  was  made  of  the  formatted  ALGOL  WRITE  statement.  This  allows 
literals  and  variables  to  appear  upon  the  same  line  of  print.  It  also 
allows  multiple  use  of  a  format  by  different  WRITE  statements.  These  are 
non-executable  statements  that  are  placed  with  the  other  such  declaritive 
statements  at  the  front  end  of  the  procedure.  A  programming  note:  it  is 


necessary  to  activate  each  format  statement  at  its  end  with  an  activation 
code  that  controls  line  feeds.  Failure  to  include  the  activation  rode 
will  have  the  effect  of  deleting  the  WRITE  statement  execution.  In  this 
sense,  the  format  statements  are  entirely  different  from  their  counter¬ 
parts  in  FORTRAN.  See  reference  [4]  for  more  details. 

.  In  addition,  a  great  many  changes  were  made  to  the  sections  of  the  processor, 
containing  the  code  for  the  individual  commands  but  these  will  not  be 
covered  here  since  this  was  already  done  in  Section  B. 
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